# Geometrically Nonlinear Static Analysis of Edge Cracked Timoshenko Beams Composed of Functionally Graded Material

Authors: Şeref Doğuşcan Akbaş
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

Geometrically nonlinear static analysis of edge cracked cantilever Timoshenko beams composed of functionally graded material (FGM) subjected to a nonfollower transversal point load at the free end of the beam is studied with large displacements and large rotations. Material properties of the beam change in the height direction according to exponential distributions. The cracked beam is modeled as an assembly of two subbeams connected through a massless elastic rotational spring. In the study, the finite element of the beam is constructed by using the total Lagrangian Timoshenko beam element approximation. The nonlinear problem is solved by using incremental displacement-based finite element method in conjunction with Newton-Raphson iteration method. The convergence study is performed for various numbers of finite elements. In the study, the effects of the location of crack, the depth of the crack, and various material distributions on the nonlinear static response of the FGM beam are investigated in detail. Also, the difference between the geometrically linear and nonlinear analysis of edge cracked FGM beam is investigated in detail.


---

## Body

## 1. Introduction

Functionally graded materials (FGMs) are special composites whose composition varies continuously as a function of position along thickness of a structure to achieve a required function. FGMs are generally made of a mixture of ceramic and metal to satisfy the demand of ultrahigh-temperature environment and to eliminate the interface problems. Typically, in an FGM, one face of a structural component is ceramic that can resist severe thermal loading and the other face is metal which has excellent structural strength. FGMs consisting of heat-resisting ceramic and fracture-resisting metal can improve the properties of thermal barrier systems because cracking and delamination, which are often observed in conventional layered composites, are reduced by proper smooth transition of material properties. The technology of FGMs was an original material fabrication technology proposed in Japan in 1984 by Sendai Group. Since the concept of FGMs has been introduced in 1980s, these new kinds of materials have been employed in many engineering application fields, such as aircrafts, space vehicles, defense industries, electronics, and biomedical sectors, to eliminate stress concentrations, to relax residual stresses, and to enhance bonding strength. Because of the wide material variations and applications of FGMs, it is important to study the responses of FGM structures to mechanical and other loadings. With the increased use of FGMs, understanding the mechanical behavior and safe performance of cracked FGM structures is very important. It is known that a crack in structure elements introduces a local flexibility, becomes more flexible and its dynamic, buckling, and static behaviors will be changed. Cracks cause local flexibility and changes in structural stiffness.In recent years, much more attention has been given to the nonlinear behavior of intact FGM beams. Rastgo et al. [1] investigated the instability of curved beams made of functionally graded material under thermal loading. Agarwal et al. [2] studied the geometrically nonlinear static and dynamic responses of functionally graded beams. Li et al. [3] investigated thermal postbuckling of functionally graded clamped-clamped Timoshenko beams subjected to transversely nonuniform temperature. Based on Kirchhoff’s assumption of straight normal line of beams and considering the effects of the axial elongation, the initial curvature and the stretching-bending coupling on the arch deformation, geometrically nonlinear governing equations of functionally graded arch subjected to mechanical, and thermal loads are derived by Song and Li [4]. Kang and Li [5] used large and small deformation theories to find the nonlinear solutions to a functionally graded material cantilever beam subjected to an end force. Kang and Li [6] analyzed the large deformation of a nonlinear cantilever functionally graded beam. Thermal postbuckling behavior of uniform slender FG beams is discussed independently using the classical Rayleigh-Ritz (RR) formulation and the versatile finite element analysis based on the von Karman nonlinear strain approximation by Anandrao et al. [7]. Kocatürk et al. [8] examined the full geometrically nonlinear static analysis of a cantilever Timoshenko beam composed of functionally graded material under a nonfollower transversal uniformly distributed load. Fallah and Aghdam [9] studied large amplitude free vibration and post-buckling analysis of functionally graded beams resting on nonlinear elastic foundation subjected to axial force by using Euler-Bernoulli beam theory and von Karman nonlinear strain approximation. Almeida et al. [10] presented geometric nonlinear analysis formulation for beams of functionally graded cross sections by means of a total Lagrangian formulation. Kocatürk and Akbaş [11] studied post-buckling analysis of functionally graded Timoshenko beam subjected to thermal loading by using the total Lagrangian Timoshenko beam element approximation. In recent years, the mechanical behavior of cracked FGM beams has been a topic of active research. Chan et al. [12] presented a displacement-based integral equation formulation for the mode 3 crack problem in a nonhomogeneous medium with a continuously differentiable shear modulus. Yıldırım [13] used the finite element method to find the solutions to the thermal shock fracture in an edge cracked functionally graded. Sridhar et al. [14] developed an effective pseudospectral finite element method for wave propagation analysis in anisotropic and inhomogeneous structures with or without vertical and horizontal cracks. Briman and Byrd [15] studied the effect of damage on free and forced vibrations of a functionally graded cantilever beam. Yang et al. [16] examined the free and forced vibration of inhomogeneous Euler-Bernoulli beams containing open edge cracks. Yang and Chen [17] discussed free vibration and buckling analysis of FGM beams with edge cracks by using Bernoulli-Euler beam theory and the rotational spring model. Free vibration and elastic buckling of beams made of FGM containing open edge cracks are studied within Timoshenko beam theory by Ke et al. [18]. Yu and Chu [19] studied the transverse vibration characteristics of a cracked FGM beam by using the p-version of finite element method. Ke et al. [20] used Ritz method to find solutions to the post-buckling behavior of FGM beams with an open edge crack based on Timoshenko beam theory and von Karman nonlinear strain approximation. Matbuly et al. [21] studied the free vibration analysis of a cracked FGM beam resting on a Winkler-Pasternak foundation by using differential quadrature method. Ferezqi et al. [22] presented an analytical investigation of the free vibrations of a cracked Timoshenko beam made up of FGM. Yan et al. [23] analyzed dynamic response of FGM beams with an open edge crack resting on an elastic foundation subjected to a transverse load moving at a constant speed. Yan et al. [24] examined the parametric instability of functionally graded beams with an open edge crack subjected to an axial pulsating excitation which is a combination of a static compressive force and a harmonic excitation force. Akbaş [25] presented static analysis of an edge cracked FGM beam resting on Winkler foundation by using finite element method. Yan et al. [26] discussed the nonlinear flexural dynamic behavior of a clamped Timoshenko beam made of FGM with an open edge crack under an axial parametric excitation which is a combination of a static compressive force and a harmonic excitation force based on Timoshenko beam theory and von Karman nonlinear strain approximation. Wei et al. [27] studied the free vibration of cracked FGM beams with axial loading, rotary inertia, and shear deformation by using an analytical method.The above review clearly indicates that nonlinear analysis of edge cracked FGM beams is very limited. Also, it is seen from the literature that the nonlinear studies of edge cracked FGM beams are investigated within von Karman nonlinear strain approximation in which full geometric nonlinearity cannot be considered. In von Karman nonlinear strain approximation, because some components of strain are neglected, satisfactory results can be obtained only for large displacements but moderate rotations. In the present study, geometrically nonlinear static analysis of edge cracked FGM Timoshenko beams under a nonfollower transversal point load at the free end of the beam is considered by using the total Lagrangian finite element method by taking into account full geometric nonlinearity. There is no restriction on the magnitudes of deflections and rotations in contradistinction to von Karman nonlinear strain approximation of the beam.In the literature studies, the effect of the parameter of FGM property distribution on the behaviour of mechanics and crack is limited in numerical results. To learn about more realistic mechanical behaviour of the FGM cracked beam, different material distributions must be used in numerical results. Hence, the effects of different material distributions on the crack are investigated in detail.The development of the formulations of general solution procedure of nonlinear problems follows the general outline of the derivation given by Zienkiewicz and Taylor [29]. The related formulations of geometrically nonlinear static analysis of edge cracked FG Timoshenko beam subjected to a nonfollower transversal point load at the free end of the beam are obtained by using the total Lagrangian finite element model. In deriving the formulations for geometrically nonlinear analysis of FGM Timoshenko, the total Lagrangian Timoshenko beam element formulations for homogeneous material given by Felippa [28] are used. There is no restriction on the magnitudes of deflections and rotations in contradistinction to von Karman nonlinear strain approximation of the beam. The cracked beam is modeled as an assembly of two subbeams connected through a massless elastic rotational spring. Material properties of the beam change in the thickness direction according to exponential distributions. The considered highly nonlinear problem is solved by using incremental displacement-based finite element method in conjunction with Newton-Raphson iteration method. In the study, the effects of the location of crack, the depth of the crack, and various material distributions on the displacements of the FGM beam are investigated in detail. Also, the difference between the geometrically linear and nonlinear analysis of edge cracked FGM beam is investigated in detail. The convergence study is performed for various numbers of finite elements.


## 2. Theory and Formulations

Consider a cantilever FGM beam of lengthL, width b, and height h, containing an edge crack of depth a located at a distance L1 from the left end, as shown in Figure 1. It is assumed that the crack is perpendicular to beam surface and always remains open. The beam is subjected to a nonfollower transversal point load in the transverse direction as seen from Figure 1.Figure 1
A cantilever FGM beam with an open edge crack subjected to a point load and cross section.In this study, Young’s modulusE(Y) varies continuously in the height direction (Y-axis) according to exponential distributions as follows:
(1)E(Y)=E0eβY,
where E0 is the Young’s modulus at the midplane (Y=0) of the beam. β is a constant characterizing the gradual variation of the material properties along thickness direction. When β=0, the material of the beam is homogeneous. According to (1) when Y=h/2, E=EB (EB is the Young’s modulus of the bottom). When Y=-h/2, E=ET (ET is the Young’s modulus of the top).


### 2.1. Total Lagrangian Finite Element Formulations of Intact FGM Timoshenko Beams

In this study, the total Langragian Timoshenko beam element is used and the related formulations are developed by using the formulations given by Kocatürk et al. [8] which was developed for FGM beam by using the formulations given by Felippa [28] for isotropic and homogeneous beam material. In the present study, finite element model of Timoshenko beam element is developed by using a two-node beam element shown in Figure 2. Each node has three degrees of freedom: two node displacements uxi and uyi, and one rotation θi about Z-axis.Figure 2
A two-nodeC0 beam element.A particle originally located atP0(X,Y) moves to P(x,y) in the current configuration, as shown in Figure 3. The projections of P0 and P along the cross sections at C0 and C upon the neutral axis are called C0(X,0) and C(xc,yc), respectively. It will be assumed that dimensions of the beam cross section do not change and that the shear distortion γ≪1 so that cosγ can be replaced by 1 (Felippa [28]):
(2)x=xc-Y(sinψ+sinγcosψ)=xc-Y[sin(ψ+γ)+(1-cosγ)sinψ]=xc-Ysinθ,y=yc+Y(cosψ-sinγsinψ)=yc+Y[cos(ψ+γ)+(1-cosγ)cosψ]=yc+Ycosθ,
where xc=X+uXC and yc=uXC. Consequently, x=X+uXC-Ysinθ and y=uYC+Ycosθ.From now on we will call uXC and uYC simply uX and uY, respectively, so that the Lagrangian representation of the motion is as follows:
(3)[xy]=[X+uX-YsinθuY+Ycosθ],
in which uX, uY and θ are functions of X only. This concludes the reduction to a one-dimensional model, as sketched in Figure 3(b). For a two-node C0 element, it is natural to express the displacements and rotation functions as linear in between the node displacements:
(4)w=[uX(X)uY(X)θ(X)]=12[1-ξ001+ξ0001-ξ001+ξ0001-ξ001+ξ][uX1uY1θ1uX2uY2θ2]=Nu,
in which ξ=(2X/L0)-1 is the isoparametric coordinate that varies from ξ=-1 at node 1 to ξ=1 at node 2.Lagrangian kinematics of theC0 beam element with X-aligned reference configuration: (a) plane beam moving as a two-dimensional body; (b) reduction of motion description to one dimension measured by coordinate X. This figure is given by Felippa [28].
(a)
(b)The Green-Lagrange strains are given as follows Felippa [28]:
(2.1)[e]=[e1e2]=[eXX2eXY]=[(1+uX′)cosθ+uY′sinθ-Yθ′-1-(1+uX′)sinθ+uY′sinθ]=[e-Yκγ],e=(1+uX′)cosθ+uY′sinθ-1;γ=-(1+uX′)sinθ+uY′sinθ;κ=θ′,
where e is the axial strain, γ is the shear strain, and κ is curvature of the beam, uX′=duX/dX, uY′=duY/dX, θ′=dθ/dX. By assuming that the material of FGM beam obeys Hooke’s law, the second Piola-Kirchhoff stresses in the beam becomes as follows
(6)s=[sXXsXY]=[s1s2]=[s10+E(Y)e1s20+G(Y)e2]=[s10s20]+[E(Y)00G(Y)][e1e2]=s0+Ee,
where E is the modulus of elasticity and G is the shear modulus and their dependence on Y coordinate is given by (1). Using constitutive equations, axial force N, shear force V, and bending moment M can be obtained as follows:
(2.1)N=∫As1dA=∫A[s10+E(Y)(e-Yκ)]dA=N0+Axxe-Bxxκ,V=∫As2dA=∫A[s20+E(Y)e2]dA=V0+Axzγ,M=∫A-Ys1dA=∫A-Y[s10+E(Y)e1]dA=M0-Bxxe+Dxxκ,
where
(8)N0=∫A0s10dA,V0=∫A0s20dA,M0=∫A0-Ys10dA,(Axx,Bxx,Dxx)=∫AE(Y)(1,Y,Y2)dA,Axz=∫AG(Y)dA,
where Axx, Bxx, Dxx, and Axz are the extensional, coupling, bending, and transverse shear rigidities, respectively.For the solution of the total Lagrangian formulations of FGM Timoshenko beam problem, small-step incremental approaches from known solutions are used. As it is known, it is possible to obtain solutions in a single increment of the external force only in the case of mild nonlinearity (and no path dependence). To obtain realistic answers, physical insight into the nature of the problem and, usually, small-step incremental approaches from known solutions are essential. Such incremental procedures are useful to reduce excessive numbers of iterations and in following the physically correct path. In the iterations, the load is divided by a suitable number according to the value of load. The loading is divided by large numbers. After completing an iteration process, the load is increased by adding load increment to the accumulated load.In this study, small-step incremental approaches from known solutions with Newton-Raphson iteration method are used in which the solution forn+1th load increment and ith iteration is obtained in the following form:
(9)duni=(KTi)-1Rn+1i,
where (KTi) is the system stiffness matrix corresponding to a tangent direction at the ith iteration, duni is the solution increment vector at the ith iteration and n+1th load increment, and (Rn+1i)S is the system residual vector at the ith iteration and n+1th load increment. This iteration procedure is continued until the difference between two successive solution vectors is less than a selected tolerance criterion in Euclidean norm given by the following:
(10)[(duni+1-duni)T(duni+1-duni)]2[(duni+1)T(duni+1)]2≤ζtol.A series of successive approximations gives by the following:(11)un+1i+1=un+1i+dun+1i=un+Δuni,
where
(12)Δuni=∑k=1idunk.The residual vectorRn+1i for a finite element is as follows;
(13)Rn+1i=f-p,
where f is the vector of external forces and p is the vector of internal forces given in Appendix.The element tangent stiffness matrix for the total Lagrangian Timoshenko plane beam element is as follows which is given by Kocatürk et al. [8] and Felippa [28]:
(14)KT=KM+KG,
where KG is the geometric stiffness matrix, and KM is the material stiffness matrix given as follows by Kocatürk et al. [8] and Felippa [28]:
(15)KM=∫L0BmTSBmdX.The explicit forms of the expressions in (14) is given in Appendix. After integration of (15), KM can be expressed as follows:
(16)KM=KMa+KMc+KMb+KMs,
where KMa is the axial stiffness matrix, KMc is the coupling stiffness matrix for FG material, KMb is the bending stiffness matrix, KMs is the shearing stiffness matrix, and the explicit forms of these expressions are given in Appendix.The geometric stiffness matrixKG, Bm, and the internal nodal force vector p remain the same as given by Felippa [28] and given in Appendix.


### 2.2. Crack Modeling

The cracked beam is modeled as an assembly of two subbeams connected through a massless elastic rotational spring shown in Figure4.Figure 4
Rotational spring model.The bending stiffness of the cracked sectionkT is related to the flexibility G by the following:
(17)kT=1G.Cracked section’s flexibilityG can be derived from Broek’s approximation [30]:
(18)(1-ν2)KI2E(a)=M22dGda,
where M is the bending moment at the cracked section, KI is the stress intensity factor (SIF) under mode I bending load and is a function of the geometry, the loading, and the material properties as well. ν indicates Poisson’s ratio which is taken to be constant since its influence on the stress intensity factors is quite limited [31]. For an FGM strip with an open edge crack under bending, the analytical solution and the expression of SIF are given Yu and Chu [19] which are obtained from the data given by Erdogan and Wu [31] through Lagrange interpolation technique:
(19)KI=6Mbh2πaF(ER,ah),
where a is the crack depth, ER is the ratio of Young’s modulus of bottom and top surfaces of the beam (EB/ET), and F is an unknown function of two independent variables. The function F is can be expressed as follows [19]:
(20)F(ER,ah)=(p4[ln(ER)]3+p5(ah)+p6(ah)2p1+p2ln(ER)+p3[ln(ER)]2+p4[ln(ER)]3+p5(ah)+p6(ah)2)×(+p9(ah)+p10(ah)2+p11(ah)31+p7ln(ER)+p8[ln(ER)]2+p9(ah)+p10(ah)2+p11(ah)3)-1,
where the coefficients p1,p2,…,p10,p11=1.1732, −0.3539, 0.0289, −0.0061, 0.6625, 3.072, −0.0014, −0.0017, 1.9917, −0.3496, and −3.0982 are given by Yu And Chu [19] that are determined by fitting equation (18) based on the least square method to the numerical values of the SIF for specific material gradients and normalized crack size given by Erdogan and Wu [31].After substituting (19) into (18) and by integrating (18), the flexibility coefficient of the crack section G is obtained as follows:
(21)G=∫0a/h72π(1-ν2)aF2(ER,a/h)E(a)h4da.The spring connects the adjacent left and right elements and couples the slopes of the two FGM beam elements at the crack location. In the massless spring model, the compatibility conditions enforce the continuities of the axial displacement, transverse deflection, axial force, and bending moment across the crack at the cracked section(X=L1), that is;
(22)uX(1)=uX(2),uY(1)=uY(2),N1=N2,M1=M2.The discontinuity in the slope is as follows:(23)kT(duY(1)dX-duY(2)dX)=kT(θ1-θ2)=M1.Based on the massless spring model, the stiffness matrix of the cracked section is as follows:(24)[K](Cr)=[1G-1G-1G1G]=[kT-kT-kTkT].The stiffness matrix of the cracked section is written according to the displacement vector:(25){q}(Cr)={θ1,θ2}T,
where θ1 and θ2 are the angles of the cracked section. In addition to the crack model and by use of usual assemblage procedure, the system stiffness matrix is as follows:
(26)[K](S)=[K](T)+[K](Cr),
where [K](T) is tangent stiffness matrix for the intact FGM Timoshenko beam which is given in (14).


## 2.1. Total Lagrangian Finite Element Formulations of Intact FGM Timoshenko Beams

In this study, the total Langragian Timoshenko beam element is used and the related formulations are developed by using the formulations given by Kocatürk et al. [8] which was developed for FGM beam by using the formulations given by Felippa [28] for isotropic and homogeneous beam material. In the present study, finite element model of Timoshenko beam element is developed by using a two-node beam element shown in Figure 2. Each node has three degrees of freedom: two node displacements uxi and uyi, and one rotation θi about Z-axis.Figure 2
A two-nodeC0 beam element.A particle originally located atP0(X,Y) moves to P(x,y) in the current configuration, as shown in Figure 3. The projections of P0 and P along the cross sections at C0 and C upon the neutral axis are called C0(X,0) and C(xc,yc), respectively. It will be assumed that dimensions of the beam cross section do not change and that the shear distortion γ≪1 so that cosγ can be replaced by 1 (Felippa [28]):
(2)x=xc-Y(sinψ+sinγcosψ)=xc-Y[sin(ψ+γ)+(1-cosγ)sinψ]=xc-Ysinθ,y=yc+Y(cosψ-sinγsinψ)=yc+Y[cos(ψ+γ)+(1-cosγ)cosψ]=yc+Ycosθ,
where xc=X+uXC and yc=uXC. Consequently, x=X+uXC-Ysinθ and y=uYC+Ycosθ.From now on we will call uXC and uYC simply uX and uY, respectively, so that the Lagrangian representation of the motion is as follows:
(3)[xy]=[X+uX-YsinθuY+Ycosθ],
in which uX, uY and θ are functions of X only. This concludes the reduction to a one-dimensional model, as sketched in Figure 3(b). For a two-node C0 element, it is natural to express the displacements and rotation functions as linear in between the node displacements:
(4)w=[uX(X)uY(X)θ(X)]=12[1-ξ001+ξ0001-ξ001+ξ0001-ξ001+ξ][uX1uY1θ1uX2uY2θ2]=Nu,
in which ξ=(2X/L0)-1 is the isoparametric coordinate that varies from ξ=-1 at node 1 to ξ=1 at node 2.Lagrangian kinematics of theC0 beam element with X-aligned reference configuration: (a) plane beam moving as a two-dimensional body; (b) reduction of motion description to one dimension measured by coordinate X. This figure is given by Felippa [28].
(a)
(b)The Green-Lagrange strains are given as follows Felippa [28]:
(2.1)[e]=[e1e2]=[eXX2eXY]=[(1+uX′)cosθ+uY′sinθ-Yθ′-1-(1+uX′)sinθ+uY′sinθ]=[e-Yκγ],e=(1+uX′)cosθ+uY′sinθ-1;γ=-(1+uX′)sinθ+uY′sinθ;κ=θ′,
where e is the axial strain, γ is the shear strain, and κ is curvature of the beam, uX′=duX/dX, uY′=duY/dX, θ′=dθ/dX. By assuming that the material of FGM beam obeys Hooke’s law, the second Piola-Kirchhoff stresses in the beam becomes as follows
(6)s=[sXXsXY]=[s1s2]=[s10+E(Y)e1s20+G(Y)e2]=[s10s20]+[E(Y)00G(Y)][e1e2]=s0+Ee,
where E is the modulus of elasticity and G is the shear modulus and their dependence on Y coordinate is given by (1). Using constitutive equations, axial force N, shear force V, and bending moment M can be obtained as follows:
(2.1)N=∫As1dA=∫A[s10+E(Y)(e-Yκ)]dA=N0+Axxe-Bxxκ,V=∫As2dA=∫A[s20+E(Y)e2]dA=V0+Axzγ,M=∫A-Ys1dA=∫A-Y[s10+E(Y)e1]dA=M0-Bxxe+Dxxκ,
where
(8)N0=∫A0s10dA,V0=∫A0s20dA,M0=∫A0-Ys10dA,(Axx,Bxx,Dxx)=∫AE(Y)(1,Y,Y2)dA,Axz=∫AG(Y)dA,
where Axx, Bxx, Dxx, and Axz are the extensional, coupling, bending, and transverse shear rigidities, respectively.For the solution of the total Lagrangian formulations of FGM Timoshenko beam problem, small-step incremental approaches from known solutions are used. As it is known, it is possible to obtain solutions in a single increment of the external force only in the case of mild nonlinearity (and no path dependence). To obtain realistic answers, physical insight into the nature of the problem and, usually, small-step incremental approaches from known solutions are essential. Such incremental procedures are useful to reduce excessive numbers of iterations and in following the physically correct path. In the iterations, the load is divided by a suitable number according to the value of load. The loading is divided by large numbers. After completing an iteration process, the load is increased by adding load increment to the accumulated load.In this study, small-step incremental approaches from known solutions with Newton-Raphson iteration method are used in which the solution forn+1th load increment and ith iteration is obtained in the following form:
(9)duni=(KTi)-1Rn+1i,
where (KTi) is the system stiffness matrix corresponding to a tangent direction at the ith iteration, duni is the solution increment vector at the ith iteration and n+1th load increment, and (Rn+1i)S is the system residual vector at the ith iteration and n+1th load increment. This iteration procedure is continued until the difference between two successive solution vectors is less than a selected tolerance criterion in Euclidean norm given by the following:
(10)[(duni+1-duni)T(duni+1-duni)]2[(duni+1)T(duni+1)]2≤ζtol.A series of successive approximations gives by the following:(11)un+1i+1=un+1i+dun+1i=un+Δuni,
where
(12)Δuni=∑k=1idunk.The residual vectorRn+1i for a finite element is as follows;
(13)Rn+1i=f-p,
where f is the vector of external forces and p is the vector of internal forces given in Appendix.The element tangent stiffness matrix for the total Lagrangian Timoshenko plane beam element is as follows which is given by Kocatürk et al. [8] and Felippa [28]:
(14)KT=KM+KG,
where KG is the geometric stiffness matrix, and KM is the material stiffness matrix given as follows by Kocatürk et al. [8] and Felippa [28]:
(15)KM=∫L0BmTSBmdX.The explicit forms of the expressions in (14) is given in Appendix. After integration of (15), KM can be expressed as follows:
(16)KM=KMa+KMc+KMb+KMs,
where KMa is the axial stiffness matrix, KMc is the coupling stiffness matrix for FG material, KMb is the bending stiffness matrix, KMs is the shearing stiffness matrix, and the explicit forms of these expressions are given in Appendix.The geometric stiffness matrixKG, Bm, and the internal nodal force vector p remain the same as given by Felippa [28] and given in Appendix.


## 2.2. Crack Modeling

The cracked beam is modeled as an assembly of two subbeams connected through a massless elastic rotational spring shown in Figure4.Figure 4
Rotational spring model.The bending stiffness of the cracked sectionkT is related to the flexibility G by the following:
(17)kT=1G.Cracked section’s flexibilityG can be derived from Broek’s approximation [30]:
(18)(1-ν2)KI2E(a)=M22dGda,
where M is the bending moment at the cracked section, KI is the stress intensity factor (SIF) under mode I bending load and is a function of the geometry, the loading, and the material properties as well. ν indicates Poisson’s ratio which is taken to be constant since its influence on the stress intensity factors is quite limited [31]. For an FGM strip with an open edge crack under bending, the analytical solution and the expression of SIF are given Yu and Chu [19] which are obtained from the data given by Erdogan and Wu [31] through Lagrange interpolation technique:
(19)KI=6Mbh2πaF(ER,ah),
where a is the crack depth, ER is the ratio of Young’s modulus of bottom and top surfaces of the beam (EB/ET), and F is an unknown function of two independent variables. The function F is can be expressed as follows [19]:
(20)F(ER,ah)=(p4[ln(ER)]3+p5(ah)+p6(ah)2p1+p2ln(ER)+p3[ln(ER)]2+p4[ln(ER)]3+p5(ah)+p6(ah)2)×(+p9(ah)+p10(ah)2+p11(ah)31+p7ln(ER)+p8[ln(ER)]2+p9(ah)+p10(ah)2+p11(ah)3)-1,
where the coefficients p1,p2,…,p10,p11=1.1732, −0.3539, 0.0289, −0.0061, 0.6625, 3.072, −0.0014, −0.0017, 1.9917, −0.3496, and −3.0982 are given by Yu And Chu [19] that are determined by fitting equation (18) based on the least square method to the numerical values of the SIF for specific material gradients and normalized crack size given by Erdogan and Wu [31].After substituting (19) into (18) and by integrating (18), the flexibility coefficient of the crack section G is obtained as follows:
(21)G=∫0a/h72π(1-ν2)aF2(ER,a/h)E(a)h4da.The spring connects the adjacent left and right elements and couples the slopes of the two FGM beam elements at the crack location. In the massless spring model, the compatibility conditions enforce the continuities of the axial displacement, transverse deflection, axial force, and bending moment across the crack at the cracked section(X=L1), that is;
(22)uX(1)=uX(2),uY(1)=uY(2),N1=N2,M1=M2.The discontinuity in the slope is as follows:(23)kT(duY(1)dX-duY(2)dX)=kT(θ1-θ2)=M1.Based on the massless spring model, the stiffness matrix of the cracked section is as follows:(24)[K](Cr)=[1G-1G-1G1G]=[kT-kT-kTkT].The stiffness matrix of the cracked section is written according to the displacement vector:(25){q}(Cr)={θ1,θ2}T,
where θ1 and θ2 are the angles of the cracked section. In addition to the crack model and by use of usual assemblage procedure, the system stiffness matrix is as follows:
(26)[K](S)=[K](T)+[K](Cr),
where [K](T) is tangent stiffness matrix for the intact FGM Timoshenko beam which is given in (14).


## 3. Numerical Results

In the numerical examples, the linear, and the nonlinear static deflections of the beams are calculated and presented in figures for different crack locations, crack depth, material distributions, and loads. To this end, with the use of the usual assembly process, the system tangent stiffness matrix and the system residual vector are obtained by using the element stiffness matrixes and element residual vectors for the total Lagrangian Timoshenko plane beam element. After that, the solution process outlined in the previous section is used for obtaining the related solutions for the total Lagrangian finite element model of the Timoshenko plane beam element. The top surface of the beam is pure aluminum with the material parameters:E=70GPa,ν=0.33. In the numerical integrations, five-point Gauss integration rule is used. Unless otherwise stated, it is assumed that the width of the beam is b=0.5m, height of the beam is h=0.5m, and length of the beam is L=5m in the numerical results.In Figure5, vertical displacement of the free end v(5,0) is calculated for various numbers of finite elements n for the point load P=200000 kN, the location of crack L1/L=0.1, the crack depth ratio a/h=0.6, and Young’s modulus ratio ER=1.5. The convergence study is calculated for the ratio L/h=10.Figure 5
Convergence analysis for displacements at the free end of the beam.It is seen from Figure5 that the vertical displacements converge perfectly after the finite element n=120. Therefore, in the numerical calculations, the number of finite elements is taken as n=120. In the nonlinear convergence analysis, the finite element analysis is not given an error for large number of the finite elements for the ratio L/h=10. Because, Timoshenko beam theory takes into account shear deformation and making it suitable for describing the behaviour of short beams in contradistinction to Euler-Bernoulli beam theory.In order to establish the accuracy of the present formulation and the computer program developed by the author, results obtained from the present study are compared with the available results in the literature. For this purpose, the nonlinear static deflections of an isotropic and intact cantilever beam(L=1000in.(25.40m), EI=180×103kip-in2(516541Nm2)) under a nonfollower point load P at the free end of the beam are compared with data presented by Fertis [32]. The horizontal deflections of the free end are found u=-4,64m, and Fertis [32] found the horizontal deflections u=-181.67inch(4,61m), and comparison of the results shows that the present results are comparatively nearby the results of Fertis [32].To further verify the present results, the normal stresses in the geometrically linear case, the nondimensional normal stresses of a cantilever FG beam under point load at the tip of the beam are calculated and compared with those of Kang and Li [33] by inserting the material properties used in this reference with the power-law exponent k. In the study of Kang and Li [33], the ratio of the elasticity modulus is taken as Et/Eb=10 (where Et and Eb are elasticity modulus of the material at the top and bottom surfaces, resp.), and the normal stresses are nondimensionalized by σxxbh2/Q0L (where Q0 is the point load). It is clearly seen that Figure 6 is the same as Figure 8(a) of Kang and Li [33].Figure 6
Dimensionless normal stress distributions along the beam thickness at fixed end for()k=0, (---)k=0.5, (·······)k=1, and (-·-·-)k=2.In Figure7, vertical displacements of the free end v(5,0) are shown for various loads P, different crack depth ratios a/h for the crack location L1/L=0.1 from the left end, the Young’s modulus ratio ER=1.5 for the geometrically linear and the nonlinear case.Load-vertical displacements curves for different the crack depth ratiosa/h, (a) linear case (b) nonlinear case.
(a)
(b)Load-vertical displacement curves for different crack locations(L1/L) from the left end, (a) linear case (b) nonlinear case.
(a)
(b)Figure7 shows that with the increase in the crack depth (a/h), the displacements increase, as expected. This is because by increasing the crack depth ratio (a/h), the beam becomes flexible. With the increase in the load, the effects of the crack on the beam increase significantly. Also, it is seen from Figure 7 that there is a significant difference between the geometrically linear case and nonlinear case for the edge cracked beam. The increase in load causes increase in difference between the displacement values of the linear and the nonlinear solutions.In Figure8, vertical displacements of the free end v(5,0) are shown for various loads P, differentiating the crack locations L1/L from the left end, the crack depth ratio a/h=0.6 for the Young’s modulus ratio ER=1.5 for the geometrically linear and the nonlinear case.Figure8 shows that when the crack locations get closer to the left end, the displacements increase. This is because the crack at the fixed end (left end) has a most severe effect in the beam. So, when the crack locations get closer to the left end, the beam gets more flexible. The crack locations get closer to the free end, the effect of crack on the FGM beam decreases. It is observed from the results that the crack locations play very important role in the static responses of the beam.In Figure9, the effect of Young’s modulus ratio ER on the displacements (v(5,0)) of edge cracked FGM beams (L1/L=0.1) is shown for different crack depth ratios a/h for the point load (P=10000 kN) for the geometrically linear and nonlinear case.The effect of Young’s modulus ratioER on the vertical displacements (v(5,0)) for various the crack depth ratio a/h. (a) Linear case (b) nonlinear case.
(a)
(b)It is seen from Figure9 that the increase in the Young’s modulus ratio ER causes decrease in the displacements. When the Young’s modulus ratio ER increases, the elasticity modulus of the beam increases according to (1). Another result of Figure 8 is that with the increase in Young’s modulus ratio ER, the differences of the crack depth ratio a/h decrease seriously. This is because an increase in the Young’s modulus ratio ER leads to increase in the elasticity modulus and the bending rigidity. As a result, the strength of the material increases. Hence, the beam becomes strength and the effect of the crack on the FGM beam decreases. It shows that with the suitable choice of ER, the negative effects of the crack can be reduced.In Figure10, the effect of Young’s modulus ratio ER on the displacements (v(5,0)) of edge cracked FGM beams (a/h=0.6) is shown for different crack locations (L1/L) from the left end for the point load (P=10000 kN) for the geometrically linear and nonlinear case.The effect of Young’s modulus ratioER on the vertical displacements (v(5,0)) for different crack locations (L1/L) from the left end. (a) Linear case (b) nonlinear case.
(a)
(b)It is observed from Figure10 that with the increase in Young’s modulus ratio ER, the differences of the crack locations L1/L decrease seriously. It is observed from the results that the distribution of the material plays an important role on the mechanical behaviour of the FGM beam. It shows that FGM is very effective for reducing the negative influence of the cracks.Figure11 displays the effect of crack depth ratio a/h on the deflected shape of the beam for the crack location L1/L=0.2 from the left end, the Young’s modulus ratio ER=1.5, and the point P=200000 kN for the geometrically linear and the nonlinear case.The effect of crack depth ratioa/h on the deflected shape of the beam. (a) Linear case (b) nonlinear case.
(a)
(b)Figure12 displays the effect of the crack location L1/L from the left end on the deflected shape of the beam for crack depth ratio a/h=0.6, the Young’s modulus ratio ER=1.5, and the point P=200000 kN for the geometrically linear and the nonlinear case.The effect of crack locationL1/L on the deflected shape of the beam. (a) Linear case (b) nonlinear case.
(a)
(b)It is seen from Figures11 and 12 that the crack depth ratio a/h and the crack location L1/L play important role on the static response of the beam. Also, it is seen from Figures 11 and 12 that there is a significant difference between the geometrically linear and nonlinear analysis for cracked FGM beam. It shows that to learn about more realistic mechanical behaviour of the FGM cracked beam, the geometrically nonlinear analysis must be considered.In Figure13, the effect of the Young’s modulus ratio ER on the deflected shape of the beam is shown for crack depth ratio a/h=0.6, the crack location L1/L=0.3 from the left end and the point P=20000 kN for the geometrically linear, and the nonlinear case.The effect of the Young’s modulus ratioER on the deflected shape of the beam. (a) Linear case (b) nonlinear case.
(a)
(b)It is obviously observed from Figure13 that the Young’s modulus ratio ER plays an important role in the deflection of the cracked FGM beam. With the increase in Young’s modulus ratio ER, the deflection of the cracked FGM beam decreases seriously. It has been mentioned before that with the suitable choice of ER, the negative effects of the crack can be reduced.


## 4. Conclusions

Geometrically nonlinear static analysis of an edge cracked cantilever FGM Timoshenko beam subjected to a nonfollower transversal point load is investigated with large displacements and large rotations. In the study, the finite element of the beam is constructed by using the total Lagrangian Timoshenko beam element approximation. Material properties of the beam change in the thickness direction according to an exponential function. The cracked beam is modeled as an assembly of two subbeams connected through a massless elastic rotational spring. The nonlinear problem is solved by using incremental displacement-based finite element method in conjunction with Newton-Raphson iteration method. The influences of the location of crack, the depth of the crack, and various material distributions on the geometrically nonlinear and linear static response of the edge cracked FGM beam are examined in detail. Convergence studies were performed. It is observed from the investigations that the crack locations and the crack depth have a great influence on the geometrically nonlinear static response of the FGM beam. The material distributions play important role in the mechanical behaviour of the cracked FGM beam. Functionally graded materials are very effective for reducing the negative influence of the cracks. It was seen from the investigations that there is a significant difference between the geometrically linear and nonlinear analysis for cracked FGM beam. To learn about more realistic mechanical behaviour of the cracked FGM beam, the geometrically nonlinear analysis must be considered.



---

*Source: 871815-2013-05-19.xml*