ANALYSIS OF THE ASSOCIATED STRESS DISTRIBUTIONS TO THE NONLINEAR FORCED VIBRATIONS OF FUNCTIONALLY GRADED MULTI-CRACKED BEAMS

Geometrically non-linear vibrations of functionally graded Euler-Bernoulli beams with multi-cracks, subjected to a harmonic distributed force, are examined in this paper using a theoretical model based on Hamilton's principle and spectral analysis. The homogenisation procedure is performed, based on the neutral surface approach, and reduces the FG beams analysis to that of an equivalent homogeneous multi-cracked beam. The so-called multidimensional Duffing equation obtained and solved using a simplified method (second formulation) previously applied to various non-linear structural vibration problems. The curvature distributions associated to the multi-cracked beam forced deflection shapes are obtained for each value of the excitation level and frequency. The parametric study performed in the case of a beam and the detailed numerical results are given in hand to demonstrate the effectiveness of the proposed procedure, and in the other hand conducted to analyse many effects such as the beam material property, the presence of crack, the vibration amplitudes and the applied harmonic force on the non-linear dynamic behaviour of FG beams.


INTRODUCTION
Structures made of functionally graded materials (FGM) have a very wide area of industrial applications, especially those involving severe environments with high temperatures or important temperature gradients, due to the properties of this type of materials characterized by a smooth and a continuous transition in the compositional profile. However, damages in the form of cracks may be observed in these materials during their operational life as a result of various external or internal influences or even during the manufacturing process. It is known that a crack in a structure produces a local change in the stiffness altering its dynamic characteristics and increasing its vulnerability. It is therefore of a crucial importance to understand and control the structural behaviour in such situations in order to ensure the safety and durability of the structural components. Structural health surveillance and vibration analysis are the most often used techniques by researchers and engineers in order to maintain the structural integrity. As beams are very currently used as structural components, the investigation of the dynamic behaviour of functionally graded cracked beams has been an active area of research in the last few years. Numerous methods, based on analytical, semi-analytical and numerical techniques, are available in the open literature, to investigate the effects of cracks on the linear and non-linear vibration characteristics in both the free and forced regimes. Sridhar et al. used a pseudo-spectral finite element method to analyse the wave propagation in functionally graded beam (FGB) and layered structures with and without vertical or horizontal cracks [1]. Yang and Chen analytically studied the influence of an open edge crack on the free vibration and elastic buckling of FGB using the Euler-Bernoulli beam theory and the rotational spring model with different boundary conditions [2]. Later, Ke et al. studied the free vibration and elastic buckling of cracked FGB with different end conditions within the Timoshenko beam theory using the same method as that of Yang and Chen [3]. Matbuly et al. used the differential quadrature method to study the free vibration of a FG cracked Timoshenko beam resting on a Winkler-Pasternak foundation [4]. Yu and Chu used the p-version of the finite element method to study the transverse vibration characteristics of cracked FGB [5]. Ferezqi et al. presented an analytical investigation of the free vibrations of a FG cracked Timoshenko beam employing a wave approach [6]. Lien et al. used [7] an alternative approach called dynamic stiffness method to establish a new form of frequency DIAGNOSTYKA, Vol. 22, No. 1 (2021) Chajdi M, Adri A, El Bikri K, Benamar R.: Analysis of the associated stress distributions to the … 102 equation of free vibration multi-cracked FG Timoshenko beams taking into account the actual position of neutral plane. Shabani and Cunedioglu [8] investigated the free vibrations of a cantilever symmetric FG sandwich Timoshenko beam with two cracks using the finite element method. Very recently, Kou and Yang developed a mesh-free boundary-domain integral equation method to analyse the free vibration behaviours of the FGB with edged cracks [9]. On the other hand, many researches have examined the case of a cracked FGB subjected to external loadings in the form of moving loads or axial forces [10][11][12][13][14][15][16]. Lately, Lin-Feng Zhu et al. [17] analysed the vibrational power flow of a cracked FG Timoshenko beam under a concentrated transverse harmonic force using the wave propagation approach to obtain the input power flow and transmitted power flow and applying the continuous wavelet transform to deal with the inverse problem of crack identification. It should be pointed out that all of the aforementioned investigations were limited to linear analysis. In terms of non-linear analysis of FGB, very few works have been devoted to non-linear free vibrations of cracked FGB. Kitipornchai et al. investigated the non-linear free vibration of FG cracked Timoshenko beam using the Ritz method and the rotational spring model [18]. Akbaş used the total Lagrangian finite element method to investigate the geometrically non-linear vibration of an edge cracked Timoshenko beam [19] and more recently Panigrahi and Pohit studied the free vibration of a cracked FGB, based on the neutral surface approach, the Ritz method and an iterative technique [20]. To the author's best knowledge, the researches regarding non-linear forced vibration responses of cracked FGB are very limited. Yan et al. studied the non-linear dynamic response of FG Timoshenko beam with an edge crack under a parametric excitation, combining a static compressive load and a harmonic excitation force [21]. Panigrahi and Pohit studied the non-linear dynamics of a FG cracked Timoshenko beam under an excitation force using the harmonic balance method in conjunction with an iterative technique [22]. The readers interested in more details on this topic can refer to [23]. The present article is an extension of our previous work [24], based on the single mode approach. The objective was particularly focused on the prediction of the nonlinear frequency response curves in the neighbourhood of the resonances considered. The applicability of this approach shows that the nonlinear effect does not appear only via the amplitude frequency dependence but also via the amplitude dependence of the cracked FGB deflection shapes. The purpose of the present contribution is the application of the non-linear forced vibration model previously reported in [25] using a multimode response approach and combined with the homogenization procedure based on the neutral surface approach previously adopted in [26,27]. This allowed an explicit calculation of multidimensional non-linear frequency response curves, the corresponding amplitude dependence deflection shapes and the associated curvatures distributions of a multi-cracked FGB subjected to a uniformly distributed harmonic force applied over the beam length. The closed-form solutions and the transfer matrix method used previously in [28,29] are employed here and the resulting frequency equation is solved iteratively by the Newton Raphson algorithm. The improved model developed in [30] is used to obtain the multi-cracked beam non-linear deflection shapes and the dependence frequency amplitudes (backbone curves) at large vibration amplitudes. A parametric study is performed to explore the effects of the gradient of material properties, the vibration amplitudes, the crack parameters and the applied harmonic force on the non-linear forced dynamic behaviour of multicracked FGB.  ( ) 0 z = corresponds to the middle surface. It should be emphasized that the Poisson's ratio  is considered constant throughout the analysis because of its weak influence on the stress intensity factor (SIF) [11]. The index 21 k E E = is introduced as a material property ratio.

Homogenisation procedure based on the neutral surface approach
In the case of FGB, the material property changes in the thickness direction are not symmetric, and consequently, the neutral plane does not coincide with the middle plane. As a result, the bendingextensional coupling effect occurs, as was observed in [18]. The distance separating the two planes can be expressed as follows [22]:

Rotational spring model
The cracked FGB is modelled as an assembly of uniform sub-segments connected by massless rotational springs as in [31] presenting the reduced local flexibility calculated using the fracture mechanics theory. It is assumed that the crack, perpendicular to the beam surface, remains always open. The bending stiffness K  of the cracked section is related to the flexibility C by [26]: With ( ) E  is the effective elastic modulus at the crack tip defined as: crack correction function of two independent variables. Through a number of numerical experiments, the function f has been expressed as a rational function [5]: Such as:

Linear vibration analysis
Based on the Euler-Bernoulli beam theory, the longitudinal displacement U and the transverse displacement W of an arbitrary point along the x and z axes are obtained as: where u and w are the displacement of a current neutral line point. The strains are given by: The stress-strain relation given by the generalised Hooke's law is: The force and bending moment resultant expressions are: 11 A , 11 B and 11 D are the extensional, coupling and bending rigidity coefficients defined by: Taking into account Eq. (3) and integrating the coefficients in Eq. (14) leads to 11 0 B = . Neglecting the axial inertia, the equations of motion for the th i sub-beam can be derived as: Integrating Eq. (15) with respect to x and using the immovable end conditions 0 u = at 0 x = and xL = leads to: The above expression for x N corresponds to that of an equivalent isotropic beam with ( ) 11 eff Assuming harmonic vibrations, the transverse displacements in the th j span can be written as: and the continuity and compatibility conditions at the th j crack location  1 K  and s are the stiffness of the transverse and rotational springs respectively. * K is the non-dimensional local rigidity due to the crack, related to the local flexibility coefficient C of the rotational spring by : Based on the matrix transfer method as well explained previously in [29], the specified conditions lead to a homogeneous linear system, which has a non-trivial solution if its determinant is set equal to zero. The resulting equation, solved iteratively using the Newton Raphson algorithm, leads to the beam natural frequencies. The corresponding mode shapes are then calculated by the usual algebraic procedure.

Non-linear vibration analysis
The Von Kàrmàn type non-linear straindisplacement relationships are: The elastic strain energy V of the beam is given as:  (32) The kinetic energy T and the strain energy of the crack c V as mentioned in [32] are given by: The potential energy associated to the force f V is given by [33]: The harmonic force is uniformly distributed and equal to d F respectively, excites the beam th i mode by: in which the ' i as are the unknown basic function contribution coefficients. Inserting Eq. (37) into Eqs. (32)(33)(34)(35) and applying Hamilton's principle gives: and introducing the following dimensionless parameters: A set of non-linear algebraic equations is then obtained:

The multi-mode approach
The last system Eq. (45) is similar to that obtained in [33]. It can be solved by the so-called second formulation [25], in order to obtain the beam non-linear forced vibration deflection shapes and frequency response curves at large vibration amplitudes. The basic idea behind this method consists on writing the contribution vector to the non-linear mode considered as      [24,34] from the single-mode approach application.

NUMERICAL RESULTS AND DISCUSSION
The numerical calculations of the present formulation have been made with a Matlab code. Through this program, it became possible to analyse the effects of the crack parameters such as the number, the location and the depth of the crack upon the vibration characteristics of cracked FGB. To test the validity and accuracy of the present analysis, following the approximated explicit method mentioned above in the case where the force is zero, the non-linear free vibration analysis is carried out and the numerical convergence is firstly performed. The first non-linear frequency ratios of a clampedclamped FGB with an edge crack computed by the classical Euler-Bernoulli beam theory are listed in Table 1 and compared to those obtained by It is observed that the trend of the present results matched appropriately with that obtained in [18]. It should be noted that the FGB under investigation in the rest of the analysis, has slenderness ratio 20% Lh= with 0.05 m h = and that all the analysis will focus on the clamped ends obtained by taking In order to validate the proposed procedure in the case of non-linear forced vibrations, the non-linear response curves, obtained by the present model in the neighbourhood of the first non-linear mode shape, is compared in Figure 2 with that obtained previously by the exact solution in [33], for the case of C-C isotropic beam (in the present work 1 k = ) excited by a harmonic concentrated force 500 c F = applied at its middle with R I S = (radius of gyration). It can be concluded that a supreme agreement can be observed. Before investigating the non-linear forced response of the cracked FGB, as described in the general formulation, it is necessary to determine the predominant mode excited in the case of a given distributed harmonic physical force, in order to justify the use of the multimode approach in its vicinity. Therefore, the percentages of the generalised distributed forces corresponding to the first five symmetrical modes, in the case of a C-C FGB with an edge crack at its center, applied over the whole beam span, corresponding to 5 k = are listed in Table 2. It can evidently be concluded that the force predominantly excites the first mode. Fig. 3. Comparison with the linear and nonlinear frequency response curves obtained by the single mode approach and the multimode approach of a C-C FGB with three edge cracks In Figure 3, a comparison is made between the linear and the non-linear frequency response curves obtained by the single-mode approach and the multimode approach of a C-C FGB with three edge cracks equitably distributed along the beam length, in the case when the FGB is excited by a harmonic distributed force 100 d F = along the whole cracked FGB, corresponding to (1,2,3) 0.1 h  = . As can be seen, considering the effect of geometric nonlinearity, the frequency response curve is increased and tends to the right showing the hardening behaviour. Furthermore, an important conclusion that can be deduced is that the two curves obtained in the non-linear analysis are very close to each other for maximum vibration amplitudes max wR up to 1.5 Table 2. Percentages of generalised distributed forces exciting the first five symmetric modes of a C-C FGB with an edge crack at its center.  and a maximum non-dimensional amplitudes up to more than once the beam thickness max 1.5 wR = , for various values of materiel properties is depicted in Figure 4. It can be clearly observed that the nonlinear curvature distribution curve in the case of the FGB having material property 0.2 k = i.e. "top surface rich ceramic" is hardly affected and more showed exactly where the crack is located (at the beam center), compared to the FGB with 5 k = i.e. "top surface rich metal".
The associated curvature distributions to the first nonlinear deflection obtained via the present model for a C-C FGB containing three edge cracks distributed along the beam length and excited by various values of a uniformly harmonic distributed force applied along the whole of the beam for both scenarios, when the FGB having material property 5 k = and 0.2 k = respectively are plotted in Figure  5   Our simulation brings also an important results listed in Table 3, about the effect of the distributed force variation on the percentage of correction introduced by taking into account the geometric nonlinearity effect in the analysis compared to the linear case. It has been found that when the intensity of the applied force increases, the percentage of correction decreases at both, near to the clamps and the middle of the cracked FGB and remains always higher in the case of a FGB with 0.2 k = compared to a FGB with 5 k = .
The effect of the vibration amplitudes on the nonlinear curvature distribution associated to the nonlinear deflection shapes, of a C-C FGB with three cracks equitably distributed along the beam length subjected to a uniformly harmonic distributed force 100 in the two scenarios of the FGB is plotted in Fig. 7. This figure clearly shows the influence of the nonlinear effect as expected, i.e. the increase in the stress distribution in accordance with the increase of the vibration amplitude, especially at the middle of the beam and in the regions of the clamps. It is also noticeable that the shape of the curvature changes and more specifically in the middle where the crack is located for the case of a FGB with 0.2 k = compared to that of a FGB with 5 k = , which implies that the FGB becomes slightly flexible when the crack starts to propagate in the totally surface rich ceramic. As can be seen from Table 4, the percentage of correction introduced by the non-linear theory compared to the linear theory increases by increasing the vibration amplitude either in the middle of the beam or in the region of the clamps, and also shows us that there is no proportionality on the stress variation with the amplitude variation contrary to the linear case because the principle of supposition applied in this domain remains invalid and remains completely fails to predict the correct behaviour of the beam response at large vibration amplitudes. Figure 8 shows the effect of crack depth on the non-linear curvature distributions of a C-C FGB with three cracks equitably distributed along the FGB length excited by a uniformly harmonic distributed force 100 d F = vibrating with significantly large amplitude max 1.5 wR = . It is well noticed that the evolution of the crack depth is rather important and compared to that of a FGB with 5 k = , the curvature decreases at the crack location when the crack depth increases, this clearly states that this phenomenon "softening" is caused by the flexibility of the beam. On the other hand, an inverted variation of curvature is obtained here at the regions of the clamps in which the curvature increases with increasing depth in the non-linear case and decreases in the linear case. It can also be deduced from the Table 5 that the crack depth effect increases the percentage correction of the curvature variation at the clamps regions and decreases it in the middle of the beam where the crack is located. The effect of a uniformly harmonic distributed force A peculiar phenomenon can be also observed in Fig. 9, showing the effect of the crack number equitably distributed along the beam length with value of crack depth 0.3 h  = on the non-linear curvature distributions corresponding to a FGB excited by a uniformly distributed force 100 d F = . The variation of crack number is obviously noticed especially for the beam having material index 0.2 k = compared to FGB with 5 k = , and it can be concluded that the crack has more effect when it is located in the middle comparing when it is located at the points considered as nodes.

CONCLUSION
The semi-analytical approach to the non-linear dynamic response of beams based on a multi-mode analysis previously developed by Benamar and coworkers has been applied in the present study to the non-linear forced vibration of functionally graded beams with multi-cracks subjected to a harmonic force for the clamped boundary conditions case. A homogenization procedure has been proposed based on the neutral surface approach leading to a simplified formulation of the non-linear problem which is reduced to that of an equivalent isotropic homogeneous multi-cracked beam without altering the accuracy of the method of analysis. A numerical solution using a simplified approximate method (the so-called second formulation) enabled the non-linear multimode frequency response function to be obtained as a solution of the multi-dimensional Duffing equation. This permitted later to obtain the non-linear forced deflection shapes and associated curvatures of multi-cracked FGB.
A comparison of the response curve of an intact isotropic beam subjected to a concentrated force with that obtained by an exact solution is carried out in order to test the validity of the proposed method. An excellent agreement had been observed.
In addition, the accuracy and convergence of the present model and of the crack modelling using the classical Euler-Bernoulli beam theory are validated by comparison with the literature results obtained by the Timoshenko beam theory in the case of a fully clamped FGB having an edge crack. An appropriate correspondence has been obtained in the tendency of the first non-linear frequency ratios according to the increase of the crack depth.
A parametric study and numerical results are also given in the form of bending moment distributions or curvatures allowing the analysis of numerous effects such as the material property, the presence of the crack and the applied harmonic distributed force on the non-linear dynamic response which can be useful in the prediction of the fatigue life of a functionally graded beam containing multi-cracks and subjected to large vibration amplitudes, providing a quantitative estimate of the percentage correction introduced by the non-linear theory compared to the linear theory.