3D energy-based finite element elasticity approach for shear postbuckling analysis of functionally graded plates on elastic foundations

Accepted Manuscript 3D energy-based finite element elasticity approach for shear postbuckling analysis of functionally graded plates on elastic foundations M. Shariyat, K. Asemi PII: DOI: Reference:

S0263-8223(16)30636-5 http://dx.doi.org/10.1016/j.compstruct.2016.05.057 COST 7470

To appear in:

Composite Structures

Received Date: Revised Date: Accepted Date:

26 February 2016 10 May 2016 17 May 2016

Please cite this article as: Shariyat, M., Asemi, K., 3D energy-based finite element elasticity approach for shear postbuckling analysis of functionally graded plates on elastic foundations, Composite Structures (2016), doi: http:// dx.doi.org/10.1016/j.compstruct.2016.05.057

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

3D energy-based finite element elasticity approach for shear postbuckling analysis of functionally graded plates on elastic foundations M. Shariyat†a , K. Asemi*b a b

Faculty of Mechanical Engineering, K.N. Toosi University of Technology, Tehran 19991-43344, Iran. Department of Mechanical Engineering, Islamic Azad University, Tehran North Branch, Tehran, Iran.

Abstract In the present paper, the shear post-buckling behavior of the FGM rectangular plates is analyzed. Furthermore, effects of the elastic foundation on the post-buckling deformation/distortion patterns and strengths are investigated. The exact energy-based 3D theory of elasticity is implemented through incorporation of the general form of Green’s strain tensor. A fully consistent 3D Hermitian element with up to the second-order mixed derivatives continuity is utilized to satisfy continuity of the stress and distortion components at the mutual nodes of the adjacent elements. The post-buckling equilibrium paths are traced using an arc-length solution scheme in conjunction with the modified Newton-Raphson method. A comprehensive sensitivity analysis is accomplished to evaluate effects of the heterogeneity index, elastic foundation stiffness, aspect ratio, and the considered two types of practical boundary conditions on the post-buckling deformation patterns, mode jumping, and equilibrium paths. Results obtained based on post-buckling observations of properly adopted reference points reveal that mode jumping and larger distortions may occur due to increasing the post-buckling loads, the elastic foundation leads to larger number of the local deformation regions, and the buckling strengths reduce remarkably when two opposite edges are clamped. Keywords: Shear post-buckling; Three-dimensional finite element; Theory of elasticity; FGM plate; Elastic foundation; Arc length method.

1 Introduction The thin-walled regions of the structures may be prone to shear buckling or post-buckling wrinkling whenever their opposite boundaries are subjected to relative shear displacements or shear forces. Plates that are used in construction of the skin panels or infra-structure of the ground, marine, space, and aerospace vehicles and propeller shafts of the marine and ground vehicles are some typical examples. The experimental observations show that buckling may be accompanied by a substantial release of the stored strain energy [1-3]. Hence, prediction and control of the buckling, wrinkling, and post-buckling deformation patterns of the dissipative elements of the structure (e.g., the buckling initiator elements designed for dissipation of the crush energy in the vehicles or trains) may be of crucial importance. Top and bottom layers of the plate may have to be fabricated from distinct materials to meet different and sometimes, opposite design criteria. In these circumstances, usage of functionally graded materials with gradual and continuous variations of the material properties is usually preferred to using the composites. †

Corresponding author, Professor, E-mail addresses: [email protected] and [email protected] Tel.: +98 9122727199; Fax: +98 21 88674748, Web page: http://wp.kntu.ac.ir/shariyat/ . * Assistant Professor, E-mail address: [email protected]

1

The available researches on shear buckling of the rectangular plates have mainly focused on the homogeneous isotropic or composite plates. Budiansky [4] studied shear buckling of the isotropic clamped rectangular plates. Smith et al. [5] presented a Rayleigh–Ritz method for local buckling analysis of rectangular unilaterally restrained plates in pure shear, using the classical plate theory (CPT). Azikov et al. [6] studied buckling behavior of composite plates subjected to combined in-plane compression and shear loads. Loughlan [7] studied shear buckling behavior of the laminated composite plates, using the finite strip procedure and the CPT. Using the CPT, Lopatin and Korbut [8] studied shear buckling of the clamped orthotropic laminated plates. Papadopoulos and Kassapoglou [9] determined buckling loads of rectangular composite plates consisting of concentric rectangular layups under shear, using the CPT. A semi-analytical approach to the buckling analysis of laminated plates subjected to various combinations of in-plane shear, compression, and tension loads was presented by Shufrin et al. for homogeneous isotropic[10] and orthotropic [11] plates, based on the CPT. Wu et al. [12] employed the extended spline collocation method to analyze linear shear buckling of the isotropic homogeneous rectangular plates based on the CPT. Explicit expressions were developed by Tarján et al. [13] for buckling analysis of long rectangular plates with constrained edges under shear and linearly varying edge loads. Chirica and Beznea [14] studied influence of the initial delaminations on buckling behavior of composite plates under shear and compression loadings. An analytical solution was developed by Yang et al. [15] for shear buckling of anisotropic rectangular plates having arbitrary edge and corner conditions. Uymaz and Aydogdu [16] presented a three-dimensional shear buckling analysis for rectangular functionally graded plates with different boundary conditions, using the linear strain energy concept. Almost all the limited shear post-buckling analyses were accomplished for the composite plates. Kosteletos [17] investigated post-buckling behavior of rectangular plates with clamped edges under inplane shear or a combination of the in-plane loads, using the first-order shear-deformation plate theory (FSDT). Mittelstedt et al. [18] investigated post-buckling behavior of rectangular composite plates under in-plane shear by a variational method. Singh and Kumar studied post-buckling of rectangular laminates under in-plane shear loads [19] or uniaxial compression combined with in-plane shear [20]. Kim et al. [21] carried out post-buckling analysis of composite plates under pure shear. Han et al. [22] studied postbuckling behavior of laminated composite plates under combinations of in-plane shear, compression, and lateral loading, using an FSDT formulation. Gupta et al. [23] studied post-buckling behavior of laminated plates containing material damages, under uniaxial, biaxial compressive and in-plane shear loadings. Post-buckling analysis of rotationally-restrained laminated plates under combined shear and compression was presented by Chen and Qiao [24], using the bending CPT. Nie et al. [25] studies post-buckling behaviors of long asymmetrical rotationally-restrained laminated composite plates subjected to shear loads, using the CPT. Buckling and post-buckling analyses of symmetric variable angle tow composite plates under combined axial compression and shear were performed by Raju et al. [26], using the CPT. Static and dynamic analyses of plates resting on elastic foundations have mainly been accomplished based on appropriate plate theories. Ferreira et al. [27] conducted static and free vibration analysis of rectangular plates resting on a two-parameter Pasternak foundations employing the FSDT and collocation with radial basis functions. Tornabene et al. presented static and dynamic analyses for laminated doublycurved shells resting on the Winkler–Pasternak elastic foundations either by using the FSDT [28] or using the interlaminar stress recovery concept in conjunction with a general higher-order formulation [29]. A hyperbolic sine shear deformation theory was used by Neves et al. [30] for the linear buckling analysis of the sandwich plates with functionally graded skins, using Carrera’s Unified Formulation and collocation with radial basis functions. Sofiyev and Kuruoglu [31] investigated torsional vibration and buckling of cylindrical shells with functionally graded (FG) coatings surrounded by two-parameter Pasternak elastic foundations. Very rare buckling analyses have been presented so far based on the general threedimensional theory of elasticity [32,33], especially for plates resting on elastic foundations. The elastic foundation changes the buckling strengths and deformation patterns of the functionally graded plates [34,35]. Shear buckling of orthotropic FGM plates surrounded by elastic foundations was investigated by Shariyat and Asemi [36], employing the 3D theory of elasticity and a novel 3D cubic B-spline element. Asemi et al. [37] employed a 3D finite element elasticity approach to investigate buckling of

2

heterogeneous functionally graded plates under biaxial compression, shear, tension-compression, and shear-compression loads. A 3D finite element elasticity formulation was developed by Asemi and Shariyat [38] for buckling of anisotropic FGM plates under uniaxial and biaxial loading conditions. Shariyat and Asemi [39] extended this approach to FGM plates subjected to non-uniform in-plane compressions. Recently, 3D post-buckling of rectangular transversely and longitudinally graded plates under compressions and on elastic foundations was investigated by Asemi and Shariyat [40] and Shariyat and Asemi [41], respectively. No research has been published on shear post-buckling of the rectangular FGM plates to date. In the present research, shear post-buckling behavior of the FGM rectangular plates is investigated using the exact 3D theory of elasticity. Moreover, effects of the elastic foundation on the post-buckling deformations and strengths are studied. A full consistent 3D Hermitian element with up to the secondorder mixed derivatives continuity is employed to satisfy continuity of the stress and distortion components at the mutual nodes of the adjacent elements. Full Green’s strain tensor is employed. Solution of the resulting highly non-linear governing equations is accomplished by using the arc-length method in conjunction with the modified Newton-Raphson procedure.

2 Description of the displacement, strain, and stress fields A schematic of the post-buckling shear deformations of the FGM plate are shown in Fig. 1 along with the geometric parameters and the employed coordinate system. Length, width, and thickness of the plate are denoted by a, b, and h, respectively. The transverse coordinate of the plate (z), is measured from the bottom surface of the plate. The elasticity modulus (E) of the FGM plate is assumed to vary in the transverse direction according to the following power law: ( E ( z ) = Em + ( Ec − Em )( z / h ) n 1) where n is the positive-definite heterogeneity index. The subscripts c and m stand for the ceramic and metal constituent materials, respectively. In addition to the edge supports, the plate is supported by a linear elastic Winkler-type foundation.

Fig. 1 Shear post-buckling deformations as well as the employed geometric parameters and coordinate system of the FGM plate. The plate is supported by edge supports and an elastic foundation not shown. The well-known commercial finite element analysis codes (e.g., ABAQUS, NASTRAN, and ANSYS) use C0-continuous Lagrangian elements and thus, cannot establish stress, twist, slope, and distortion continuity at the mutual nodes and edges of the adjacent elements. Therefore, their results may be used

3

only when a huge number of elements is used. To avoid this shortcoming, the plate is discretized by using a full consistent 3D Hermitian brick element with up to the second-order mixed derivatives continuity. The natural coordinates −1 ≤ ξ ≤ 1 , −1 ≤ η ≤ 1 , and −1 ≤ ζ ≤ 1 are defined along the x, y, and z coordinates, respectively, to present the shape functions of the element. The following relations may be presented to relate the global derivatives to the local natural ones: ∆x−1 0 0   ∂ x  ∂ ξ           ∆y−1 0  (2) ∂ y  = 2  0 ∂ η        − 1 ∂   0 0 ∆z   ∂ z       ζ  where ∆x , ∆y , and ∆z are dimensions of the element. The displacement components u, v, and w are defined along the x, y, and z directions, respectively. The element may be mapped form the geometric space to a cubic Hermitian element with side length of 2 (in the natural coordinates) and nodal degrees of freedom that are show in Fig. 2. Variations of the displacement components within the element may be related to the nodal quantities through the following equation: u     Ω(ξ , η , ζ ) =  v  = Λ(ξ , η , ζ )Ω(e ) (3)    w  where Λ is the shape functions matrix and the nodal matrix is: Ω ( e )T = U ( e ) , V ( e ) , W ( e )

(4)

and  (e ) U = (U ,U ,ξ , U ,η ,U ,ζ , U ,ξη ,U ,ξζ , U ,ηζ )1 ,…, (U , U ,ξ ,U ,η ,U ,ζ , U ,ξη ,U ,ξζ ,U ,ηζ )8   (e ) (5) V = (V ,V,ξ ,V,η ,V,ζ , V,ξη , V,ξζ ,V,ηζ )1 ,…,(V ,V, ξ ,V,η ,V,ζ , V,ξη ,V,ξζ ,V,ηζ )8   ( e)  W = (W ,W,ξ ,W,η ,W,ζ , W,ξη ,W,ξζ ,W, ηζ )1 ,…, (W ,W,ξ ,W,η ,W,ζ , W,ξη ,W,ξζ ,W,ηζ )8  where the comma symbol denotes a partial derivative with respect to the indicated coordinate, respectively.

Fig. 2 Mapping of the original geometric element to the full consistent Hermitian cube element with 21 degrees of freedom per node. Each node of the Hermitian element has 21 nodal degrees of freedom, as shown in Fig. 2. Presence of the mixed second-order derivative degrees of freedom guarantees path-independence of the Galerkin integrals and normal derivatives of the displacement quantities. Denoting the nodes, in each direction by 1

4

and 2 if value of the mentioned coordinate is -1 and 1, respectively, one may define the shape functions matrix of Eq. (3) as: Ψ (ξ , η , ζ )  0 0    Λ(ξ , η , ζ ) =  0 Ψ( ξ , η , ζ ) 0 (6)     0 0 Ψ( ξ , η , ζ )   where the Hermitian shape functions matrix ( Ψ1×56 ) may be defined using the idea of separation of variables, due to the tabular arrangement of the nodal points, as: Ψ(ξ , η , ζ ) = N1 (ξ , η , ζ ),…, N 8 (ξ , η , ζ ) N i (ξ , η , ζ ) = Hi 0 (ξ ) Hi0 (η )Hi 0 (ζ ) Hi1 (ξ ) Hi0 (η ) Hi 0 (ζ ) Hi 0 (ξ ) Hi1 (η ) Hi 0 (ζ ) Hi 0 (ξ )Hi0 (η ) Hi1 (ζ ) Hi1 (ξ ) Hi1 (η )Hi 0 (ζ )

Hi1 (ξ ) Hi0 (η )Hi1 (ζ )

where, e.g., Hi 0 (ξ ) = 14 (ξ + ξi )2 (2 − ξξi )

(7)

0

Hi (ξ ) Hi1 (η )Hi1 (ζ )

(8)

Hi1 (ξ ) = 14 ξi (ξ + ξi ) 2 (ξξi −1)

where ξi , is value of the natural coordinate ξ associated with the i-th node. In thin, displacement-based elements (e.g., when using plate theories), full integration of the stiffness matrices may lead to shear and membrane locking or over-stiff behavior, especially when highly constrained boundaries or coarse meshes are used [42]. These problems are owing to the low number of degrees of the element; so that sides of the element cannot deform in harmony with the applied loads. Therefore, the linear (e.g., the 4-node tetrahedron), bi-linear (e.g., the 4-node quadrilateral) and tri-linear (8-node brick) elements are obviously too stiff and susceptible to locking. Therefore, higher forces/moments are required to obtain a certain displacement pattern than for the quadratic elements. Appearance of tri-linear terms in the shape functions of the 8-node brick elements leads to volumetric (incompressibility) locking [43]. Some authors have varied factors of the strain energy [44] or employed reduced/selective integration techniques as a remedy [42,43] in large deformation problems in finite elasticity. Yang et al. [44] presented a family of one-dimensional finite elements for the analysis of wrinkling in stiff thin films resting on a thick elastic substrate. While Euler–Bernoulli’s theory was used for the thin film, higher-order approximations were used for the substrate. They used a reduced form of Hooke’s law for the thin film to avoid the Poisson locking due to using the 1D Euler–Bernoulli theory. Hazrati and Shariyat [45] have proven that analytically determination of the FEM integrals may be considered as another solution for the numerical or shear locking. Quadratic shape element (e.g., in quad elements with 8 or 9 nodes) have additional degrees of freedom for the sides that allow them to curve and represent the actual deformation of the element much better [46] and thus, manifest no shear/membrane/bending locking or mechanism formation [42]. On the other hand, shear locking generally occurs when not enough elements exist across the thickness (e.g., when shell element are used) or when mid-side nodes are not used. No locking may be imagined for the present Hermitian element due to its significant flexibility, and capability to represent/trace the actual deformations owing to containing 168 degrees of freedom. Moreover, due to using the 3D elasticity, the thickness can be divided into several elements, in contrast to cases wherein plate theories (all types of the plate theories, e.g., equivalent single-layer, layerwise, or zigzag theories) are employed. To accurately trace events of the post-buckling region wherein the deformations and rotations are large, using the complete form of Green’s strain tensor instead of the von Karman assumptions is preferred. The full Green’s strain tensor for curved coordinates may be written as: 1 1 (9) [ ε ] = [∇Ω + (∇Ω)T + ∇Ω ⊗ (∇Ω)T ] = [Ωi , j ei ⊗ e j + Ωj ,i e j ⊗ ei + Ωα,iΩα, j ei ⊗ e j ] 2 2

5

where ⊗ is a dyadic multiplication notation and α is a dummy index. In the rectangular Cartesian coordinates, the tensor quantities are equivalent to their physical components; so that, Green’s straindisplacement relations become:   ε x   u, x   12 (u,2x + v,2x + w,2x )              1   εy  v, y    (u,2y + v,2y + w,2y )      2             2 2 2   1 ε w     z  ,z ( u + v + w )      , z , z , z (γ ij = 2εij ) ε= = (10)  + 2  ,      γ u + v xy  ,y ,x  u, x u, y + v, x v, y + w, x w, y               γ xz  w, x + u, z        u, x u, z + v, x v, z + w, x w, z               γ   v + w    yz , z , y u u + v v + w w          , y ,z , y ,z   , y ,z   where the ε and γ symbols denote tensile and shear strains, respectively. Therefore, the strain matrix may be decomposed into linear and non-linear components: (11) ε = ε L + ε NL = ∂ Ω + ε NL = B Ω(e) + ε NL where ∂   ξ 0 0   ∆x      ∂η  0 0    ∂ x 0 0 ∆y        0 ∂y 0  ∂ ζ   0   0   0 0 ∂ z  ∆z   ∂ =  = 2 , B = ∂Λ (12) ∂    η ∂ξ  ∂ y ∂ x 0  0   ∂   ∆y ∆x   z 0 ∂x      0 ∂ ∂ ∂ ∂  ζ  z y  ξ  0    ∆z ∆x    ∂ζ ∂η    0   ∆z ∆y  Therefore, based on Eqs. (6) and (12):  Ψ 0 0   ,ξ   0 κ1Ψ,η 0    0 κ2 Ψ,ζ  ∆x ∆x 2  0 B = , κ1 = , κ2 = (13)   Ψ,ξ 0  ∆x  κ1Ψ,η ∆y ∆z κ Ψ 0 Ψ,ξ   2 ,ζ   κ2 Ψ,ζ κ1Ψ,η   0 The non-linear components of the strain matrix can be expressed in terms of the nodal values, as follows:

6

( e )T T    ( e )T   Ω Λ, ξ Λ, ξ   Ω Πξξ   2 ( e )T T   ( e )T   κ1 Ω Λ,η Λ,η  Ω Πηη   2 ( e )T T   ( e )T  2  κ2 Ω Λ,ζ Λ,ζ  ( e)  Ω Πζζ  (e ) nl Ω = Ξ(Ω( e) )Ω(e ) ε = Ω = ( e )T Ω Π  (∆x) 2  2κ1Ω(e)T ΛT,ξ Λ, η  ξη      ( e )T  ( e )T T  2κ2 Ω Λ,ξ Λ,ζ   Ω Πξζ     ( e )T   2κ1κ2Ω(e )T ΛT,η Λ,ζ  Ω Πηζ      Since εij = ε ji , it can easily be verified that: ( e )T T    ( e )T T   Ω Λ, ξ Λ, ξ   Ω Πξξ   2 ( e )T T   ( e )T T   κ1 Ω Λ, η Λ,η  Ω Πηη   2 ( e )T T   ( e )T T     Ω Πζζ  ( e ) κ Ω Λ Λ 2 2 , ζ , ζ   Ω(e ) =  Ω εnl = 2  ( e )T T ( e )T T    2 κ Ω Λ Λ Ω Π ∆ x ( )  1 ,η ,ξ  ξη      ( e )T T  ( e )T T  2κ2 Ω Λ,ζ Λ,ξ  Ω Πξζ     ( e )T T  ( e ) T T  2κ1κ2Ω Λ,ζ Λ,η  Ω Πηζ      So that, Π is always a symmetric matrix: Πij = ΠijT

(14)

(15)

(16) Indeed, this property is due to the fact that the strain tensor is always a symmetric tensor (in contrast to the stress tensor that may be asymmetric in presence of the couple stresses or body couples). Eq. (11) may be rewritten as: (17) ε = ∂ Ω + ε NL = [ B + Ξ(Ω ( e ) )]Ω ( e ) The stress and strain matrices may be related through Hooke’s generalized law: (18) σ = Q ε = Q [ B + Ξ(Ω ( e ) )]Ω ( e ) where Q is the stiffness matrix and the stress vector is defined as:

σT = σ x

σy

σz

τ xy

τ xz

(19)

τ yz

3 Governing equations of the shear post-buckling behavior The governing equations of the shear post-buckling behavior are derived based on the principle of minimum total potential energy: δ Π = δU − δW = 0 (20) where U and W are the strain energy and work of the externally applied loads, respectively, and defined as: δU =

∫

δεT σdV

(21)

V

δW =



∫ (s δΩ) Γ

T

T

 τ  dΓ − 

∫

A

K w wδ wd A

(22)



in which, τ is the shear traction acting on the edge of the plate, s is the tangent unit vector of the boundary, K w is coefficient of the Winkler-type elastic foundation, and V, A , and Γ are the volume, area, and edge surface, respectively. According to Eqs. (11), (14), (16), and noting that δεij = δε ji , one may deduce that:

7

  ( e )T   ( e )T      δΩ Πξξ   Ω Πξξ     ( e )T   ( e )T     δΩ Πηη  Ω Πηη      ( e )T   ( e )T      δΩ Πζζ  (e )    Ω Πζζ   (e )  (e )  Ω  =  B + 2   δΩ = [ B + 2Ξ(Ω(e ) )]δΩ(e ) δε =  B δΩ + 2  (e )T (23) ( e )T  δΩ Π         ξη    Ω Πξη       ( e )T   ( e )T      δ Ω Π Ω Π      ξζ ξζ         ( e )T ( e )T         δ Ω Π Ω Π ηζ  ηζ        On the other hand, Eq. (23) implies that: (24) δΞ(Ω( e ) )Ω ( e ) = Ξ(Ω( e ) )δΩ ( e ) Therefore, according to Eqs. (18) and (23), Eq. (21) may be rewritten as:   δU = δΩ( e )T  {B TQB + B TQΞ(Ω( e ) ) + 2ΞT (Ω( e ) )Q[ B + Ξ(Ω( e ) )]} dV  Ω( e ) (25) V   where: ∆x∆y∆z dV = dxdydz = d ξ d ηd ζ (26) 8 On the other hand, according to Eq. (22): 0 0 0   ( e )T T  (e)T T  δW = (δΩ Λ sτ ) d Γ − K w δΩ Λ 0 0 0 ΛΩ( e ) d A (27)   A 0 0 1 Γ   Because:  0 0 0   T (e) T (e) ( e )T T  δ w w = ([0 0 1]ΛδΩ ) ([0 0 1]ΛΩ ) = δΩ Λ 0 0 0 ΛΩ( e ) (28)    0 0 1   Substituting Eqs. (25) and (27) into Eq. (20), and recalling that generally, δΩ(e ) ≠ 0 at the interior regions of the plate, yields the following finite element governing equations for the post-buckling deformations:    0 0 0     ( ) ( ) ( ) T T e T e e T  {B QB + B QΞ(Ω ) + 2Ξ (Ω )Q[ B + Ξ(Ω )]} dV + K Λ 0 0 0 Λd A  Ω( e ) w  V   A    0 0 1 (29)      = (ΛT sτ ) d Γ

∫

∫

∫

∫

∫

∫ Γ

Eq. (29) can be represented by the following compact form: K L + K NL (Ω( e ) ) Ω( e ) = F   where Kl and K NL are the linear and non-linear stiffness matrices: 0 0 0   L T T  K = B QB dV + K w Λ 0 0 0 Λd A   V A 0 0 1  

∫

K NL (Ω( e ) ) =

∫

∫ {B V

T

QΞ(Ω( e ) ) + 2ΞT (Ω( e ) )Q [ B + Ξ(Ω( e ) )]} dV

4 The combined arc-length and modified Newton-Raphson solution procedure

8

(30)

(31)

(32)

Since the post-buckling governing equations are non-linear, they be solved through incrementally increasing the applied loads and determination of the resulting deflection [47,48]. Due to the complicated curvature variations of the post-buckling curve, such as, existence of the snap-through and snap-back segments, the common iterative methods may fail in tracing the load-displacement equilibrium paths and passing the limit points; because the load level is held constant whilst iterating to convergence. Hence, a load scaling factor called load parameter has to be defined to adjust the next load level associated with the post-buckling equilibrium. This task may be fulfilled through using the arc-length technique [49-52]. Some researchers have employed the arc-length technique for tracing the post-buckling behaviors [5357]. In the present research, a combination of the Crisfield-Ramm arc-length and modified NewtonRaphson techniques is used to relate the solution associated with the i+1 load step to that of the ith loading step. Since resultants of the internal nodal forces have to be equal to the external ones; so that: ( i +1 1 i +1 R = i +F F = Fi + i+1 ∆F , 33) According to the tangent stiffness matrix of the previous loading step, one may write: ( i +1 1 i +1 i i +1 ∆F = i +K δΩ ( e ) ≅ K δΩ ( e ) T T 34) On the other hand, the tangent stiffness matrix may be written based on Eq. (30) as: i    ∂K NL  (e )   ∂i F ∂i R ( (e )  i L NL   Ω  K T = i (e ) = i (e) =  K + K ( Ω ) +    ( e )   35)     ∂Ω ∂Ω  ∂Ω      where according to Eq. (32): ∂K NL ( =  B TQΞ(Ω( e) ) + 2[Ξ(Ω( e ) )]T Q ( B + 2[Ξ(Ω( e) )]) dV (e )   36) V ∂Ω Substituting Eqs. (34) into Eq. (33) gives ( i 1 K T i+1δΩ ( e) = i+R − Fi 37) The nodal forces can be assumed to be multiplier of a reference force F 0 ; so that: ( i +1 i +1 λ = iλ + i ∆λ Fi = iλF 0 , F = i+1λF 0 , 38) where λ is the so called load parameter. Denoting the iteration counter by j, Eq. (37) may be rewritten as follows for the jth iteration, based on Eq. (38): ( i 1 K T i+j1δΩ( e ) = i+R − iλF 0 − i+1j ∆λF 0 = ∆P ( iλ) − i+1j ∆λF 0 39) Therefore, i+j1δΩ( e ) may be determined based on an iterative solution. According to Eq. (37), Eq. (39) may

∫

be solved for i +1 ( e) j δΩ

i +1 ( e) j δΩ

as:

( 40)

= i K T −1∆P ( iλ ) − i +1j ∆λΩT(e)

where

( 41) The total incremental displacements and load parameter of the current (i+1) step may be determined as follows: i +1 δΩ(e ) = iδΩ(e ) + i +j1δΩ( e) ( ΩT( e ) = i K T −1F 0

42)

i +1

λ = iλ + i+1j ∆λ

Variations of the load factor ( ∆λ ) and consequently,

i +1 ( e) j δΩ

, are governed by the jth arc-length

increment that usually satisfies the following constrain proposed by Crisfield [53]:

9

( e )T i+1 (e ) i +1 j δΩ j δΩ

( 43)

= ∆li2+1

where the arc-length ∆li+1 may be updated as follows for each loading step:

( 44) where  d is the desired number of iterations for the current load step and  i is the number of iterations corresponding to convergence of the ith loading step. ∆li+1 = ∆li  d /  i

5 Results and conclusions 5.1 General specifications of the considered plates As mentioned before, the considered plate is a transversely graded one whose elastic modulus varies according to Eq. (1). Before solving the governing equations, the boundary conditions have to be imposed. The traction edge conditions have already been incorporated in the functional of the total potential energy of the system, due to incorporating their relevant works. The condensation technique outlined by Asemi and Shariyat [38] is employed here to replace the nodal second-order mixed derivatives of the displacement parameters by their equivalents in terms of the lower-order nodal quantities. The following sets of boundary conditions are adopted to extract the results: i) All edges are simply supported: x = 0, a : W ,W, η = 0; y = 0, b : W ,W,ξ = 0 (45) However, to avoid the rigid body motions, the mid-point of the plate has to be set as a reference for the in-plane displacements; so that, the other in-plane displacements will be measured relative to that point. ii) Edges parallel to the y-axis (x=0, a) are clamped and the remaining edges are free: x = 0, a : U ,V ,W ,W,η ,W,ζ = 0 (46) The base geometric and material specifications of the plate are chosen as: a / b = 1, h / a = 0.1, Ec = 151GPa, Em = 70 GPa, ν = 0.3 To present more general results, the following dimensionless load and lateral deflection are defined: w 100 P PN = , Wn = (47) Ec h   Where P = τ is magnitude of the shear stress and τ = -P (i.e., direction of the shear stresses are opposite to that of Fig. 1). 5.2 Verification of the shear post-buckling results Shear post-buckling of the FGM plate has not been investigated so far. For this reason, present results are verified by the 3D elasticity results of the ANSYS Workbench commercial finite element analysis code. However, since ANSYS cannot model the FGM plate directly, the plate has to be divided into many isotropic and homogeneous layers and consequently, the results contain some errors that depend on number of the discretization layers. In this regard, post-buckling curve of a movable simply supported FGM square plate subjected to shear loads is extracted based on the present formulation and compared with that of ANSYS Workbench. 3D 8-node HEX8 brick elements are used in ANSYS to discretize the solution domain. The thickness is divided into 10 isotropic homogenous layers with identical thicknesses and a 20×20×10 mesh is employed (Fig. 3). Results of ANSYS are plotted in Fig. 4 versus different mesh sizes. Further transverse discretization of the thickness may lead to numerical locking. While results of ANSYS exhibit a degree of convergence in the pre-buckling region, the case is somewhat different in the post-buckling region, owing the nature of the solution procedures employed by ANSYS computer code.

10

Since order of the proposed Hermitian element is significantly higher than that of ANSYS element, a 15×15×5 mesh is adopted (based on a convergence analysis whose results are not reported here), for extraction of present results.

Fig. 3 Discretization of the FGM plate in ANSYS, after dividing the thickness into 10 sub-layers.

15

PN

10

5 Present ANSYS ANSYS ANSYS ANSYS

0 0

0.1

(x=a/2, y=b/2, z=h/2) WORKBENCH 20*20*10 elements (x=a/2, y=b/2, z=h/2) WORKBENCH 20*20*8 elements (x=a/2, y=b/2, z=h/2) WORKBENCH 20*20*6 elements (x=a/2, y=b/2, z=h/2) WORKBENCH 20*20*4 elements (x=a/2, y=b/2, z=h/2)

0.2

0.3

0.4

0.5

WN Fig. 4 Results of ANSYS computer code for different mesh sizes.

11

Fig. 5 A comparison between present and ANSYS results for shear post-buckling of the movable simply supported FGM square plate (n=1) Present post-buckling results are compared in Fig. 5 with those of ANSYS, for n=1 based on the dimensionless load and lateral deflection defined in Eq. (46). A good concordance can be seen between these results. However, results of ANSYS have experienced a premature divergence. This issue may have occurred due to using C0-continuous elements (non-consistent stress field) and von Karman straindisplacement relations in ANSYS, in addition to the approximation employed in modeling the FGM material in ANSYS. Unfortunately, ANSYS Workbench computer code does not support the arc-length method in solving the nonlinear problems. Therefore, ANSYS Workbench has convergence difficulties in the deep post-buckling regions. The Full Newton-Raphson method with active Line Search option was employed to verify our results. To ensure that present results are accurate, the iterations within each load increment were terminated when the maximum relative differences between results of successive iterations became less than 0.01%. Post-buckling deformations of the plate predicted by ANSYS are illustrated in Fig. 6. Since the plate is under pure shear, the compressive principal stress acts in a direction that makes an angle of 45 degrees with the edges. For this reason, the deformation wave has formed in that direction, as Fig. 6 shows.

12

Fig. 6 Transverse post-buckling deformations pattern of the FGM plate according to ANSYS (n=1, PN = 12.12 ). 5.3 Influence of various parameters on shear post-buckling of an FGM plate without elastic foundation In the present section, a sensitivity analysis is accomplished for an FGM plate without elastic foundation. Specifications of the plate are as those introduced in Sec. 5.1. The results are extracted for a plate with either all-edges simply supported or two (opposite) clamped edges. Two reference points are chosen to trace the shear post-buckling behavior of the plate: (x=a/2, y=b/2, z=h/2) and (x=a/4, y=b/2, z=h/2). The relevant dimensionless results are plotted in Figs. 7 and 8, respectively, for various heterogeneity indices. As may be expected, due to the materials asymmetry and consequently, the extension-bending coupling, buckling strength of the FGM plates are lower than that of a plate made of a homogeneous material (i.e., the pure-metallic or pure-ceramic plates), as Figs. 7 and 8 show. Only the load-deformation of the homogeneous (i.e., pure metallic) plate includes a sharp buckling corner (bifurcation point). However, the post-buckling stiffness and strength grow with the volume fraction of the ceramic phase (i.e., by decreasing the heterogeneity indices). Since the post-buckling curve of the he pure-metallic plate is almost flat, this plate may experience an abrupt reduction in the stiffness after buckling. Since the trends observed by the two reference points agree, it can be deduced the both references points are located on the same post-buckling deformation wave. This evidence may be checked through Fig. 9 that illustrates the post-buckling transverse deformation pattern of the mid-plane of the FGM plate for n = 1 and PN = 13.91 . Since the considered plate is square, only one post-buckling deformation wave has been formed. The deformation pattern confirms that proper reference points are chosen for tracing the post-buckling behavior [40].

13

15

PN

10

5

0 0

Pure Metal n=1 n=3 n=5

0.1

0.2

0.3 WN

0.4

0.5

0.6

Fig. 7 Shear pre- and post-buckling behaviors of the simply supported square FGM plate for various heterogeneity indices, according to deflections measured at (x=a/2, y=b/2, z=h/2).

15

PN

10

5 Pure Metal n=1 n=3 n=5

0 0

0.05

0.1

0.15 WN

0.2

0.25

0.3

Fig. 8 Shear pre- and post-buckling behaviors of the simply supported square FGM plate for various heterogeneity indices, according to deflections measured at (x=a/4, y=b/2, z=h/2).

14

1 0.3 0.8 0.25 0.6 y (m)

0.2 0.15

0.4

0.1 0.2

0 0

0.05 0.2

0.4 0.6 x (m)

0.8

1

0

Fig. 9 Post-buckling transverse deformations of the mid-plane of the simply supported square FGM plate ( n = 1, PN = 13.91 ). As a second analysis, effects of the aspect ratio and choice of the proper reference points on the observed post-buckling results are investigated. In this regard, a rectangular plate with a/b=2 aspect ratio and the material properties mentioned in Sec. 5.1, is considered. Two reference points are adopted in this regard: (x= a/2, y=b/2, z=h/2) and (x=a/4, y=b/2, z=h/2) and the relevant load-deflection curves are shown in Figs. 10 and 11, respectively, for different heterogeneity indices. Although the post-buckling characteristics (e.g., the buckling loads and post-buckling stiffnesses) observed by the two reference points are identical, comparing Figs. 10 and 11 shows that these points are located on different and opposite deformation waves. Figs. 12 to 14 that are plotted for various dimensionless loads (n=1), confirm this evidence. Since buckling mechanism of the rectangular plate is quite different from that of the square plate, in contrast to the previous example (square plate), buckling strengths of the FGM plate are higher than that of the pure-metallic plate as well. Comparing Figs. 12 to 15 reveals that at loads that are slightly higher than the buckling load, the plate buckles into two deformation waves that make approximately 45 degrees angles with respect to the y-axis (i.e., the first buckling mode is associated with two buckling waves). But as the shear load increases, the wave number increases to 3, the deformation waves rotate, and the difference between angles of the deformation waves increases. In the latter situations, the mean angle of the deformation waves with respect to the y-axis is about 30 degrees. Therefore, a mode jumping has occurred. Comparing Figs. 13 and 14 implies that peaks of the first and last deformation waves has moved toward the top and bottom edges, respectively (so that, larger distortions have taken place) in higher shear loads.

15

12 10

PN

8 6 4 Pure Metal n=1 n=3 n=5

2 0

-1

-0.5

0 WN

0.5

1

Fig. 10 Shear post-buckling behavior of the simply supported rectangular (a/b=2) FGM plate for various heterogeneity indices, according to deflections measured at (x= a/2, y=b/2, z=h/2).

12 10

PN

8 6 4 Pure Metal n=1 n=3 n=5

2 0 -0.5

0 WN

0.5

Fig. 11 Shear post-buckling behavior of the simply supported rectangular (a/b=2) FGM plate for various heterogeneity indices, according to deflections measured at (x=a/4, y=b/2, z=h/2).

16

1 0.035 0.8

0.03 0.025

y (m)

0.6 0.02 0.4

0.015 0.01

0.2 0.005 0 0

0.5

1 x (m)

1.5

2

0

Fig. 12 Post-buckling transverse deformations of the mid-plane of the simply supported rectangular (a/b=2) FGM (n=1) plate, for PN = 7.95 .

1 0.1 0.05

0.8

0 -0.05

y (m)

0.6

-0.1 0.4

-0.15 -0.2

0.2

-0.25 -0.3

0 0

0.5

1 x (m)

1.5

2

Fig. 13 Post-buckling transverse deformations of the mid-plane of the simply supported rectangular (a/b=2) FGM (n=1) plate, for PN = 8.74 .

17

1

0.6 0.4

0.8 0.2 0

y (m)

0.6

-0.2 0.4

-0.4 -0.6

0.2

-0.8 0 0

0.5

1 x (m)

1.5

2

Fig. 14 Post-buckling transverse deformations of the mid-plane of the simply supported rectangular (a/b=2) FGM (n=1) plate, for PN = 11.72 . As a next step, effects of the plate thickness on the shear post-buckling behavior is investigated. In this regard, the same rectangular FGM plate is considered but the thickness is halved. The same reference points are chosen and the relevant shear post-buckling curves are plotted in Figs. 15 and 16, respectively. The resulting buckling loads are remarkably smaller than those of the previous case; but the reduction is neither proportional to the thickness nor its square. Although the post-buckling stiffness is smaller than that of the previous case, the stiffness recovery is apparent in the post-buckling regions of Figs. 15 and 16. The post-buckling deformation patterns of the mid-plane of the FGM plate are shown in Figs. 17 to 19, for various magnitudes of the dimensionless shear load. The conclusions drawn from Figs. 12 to 14 are confirmed once again; so that the mode jumping and distortion increase with load may be seen in Figs. 17 to 19 as well.

18

7 Pure Metal n=1 n=3 n=5

6 5

PN

4 3 2 1 0 -1.5

-1

-0.5

0 WN

0.5

1

1.5

Fig. 15 Shear post-buckling behavior of the thinner (h/a=0.025) simply supported rectangular (a/b=2) FGM plate for various heterogeneity indices, according to deflections measured at (x= a/2, y=b/2, z=h/2).

7 6 5

PN

4 3 2

Pure Metal n=1 n=3 n=5

1 0

-1

-0.5

0 WN

0.5

1

Fig. 16 Shear post-buckling behavior of the thinner (h/a=0.025) simply supported rectangular (a/b=2) FGM plate for various heterogeneity indices, according to deflections measured at (x=a/4, y=b/2, z=h/2).

19

x 10 7

1

-3

6 0.8 5 4

y (m)

0.6

3 0.4 2 1

0.2

0 0 0

0.5

1 x (m)

1.5

2

Fig. 17 Post-buckling transverse deformations of the mid-plane of the simply supported rectangular FGM (n=1) plate, for PN = 2.32 (a/b=2, h/a=0.025, n=1).

1

0.3 0.2

0.8

0.1 0

y (m)

0.6

-0.1 -0.2

0.4

-0.3 0.2

-0.4 -0.5

0 0

0.5

1 x (m)

1.5

2

Fig. 18 Post-buckling transverse deformations of the mid-plane of the simply supported rectangular FGM plate, for PN = 2.98 (a/b=2, h/a=0.025, n=1).

20

1

1.5 1

0.8

0.5

y (m)

0.6

0 0.4 -0.5 0.2

0 0

-1

0.5

1 x (m)

1.5

2

Fig. 19 Post-buckling transverse deformations of the mid-plane of the simply supported rectangular FGM plate, for PN = 6.62 (a/b=2, h/a=0.025, n=1). Finally, effects of the edge conditions are investigated through studying post-buckling behavior of an FGM plate whose x=0 and a edges are clamped. The geometric and material specifications are as those mentioned in Sec. 5.1. Once again, the (x= a/2, y=b/2, z=h/2) and (x=a/4, y=b, z=h/2) references points are chosen for measurement of the lateral deflection. The relevant shear post-buckling curves are illustrated in Figs. 20 and 21, respectively, for different heterogeneity indices. While these figures confirm each other, they show a higher post-buckling stiffnesses in comparison to those of Figs. 7 and 8, respectively. The post-buckling deflection pattern of the FGM plate (n=1) is plotted in Fig. 22, for a specific non-dimensional shear load. Fig. 22 shows that in contrast to the simply supported plate (Fig. 9), even though the plate is square, the first buckling mode corresponds to two distinct transverse deformation waves. In contrast to the simply supported plate wherein the deformation wave was formed at the middle of the plate, these deformation waves are located in the neighborhood of the corners of the free edges. However, in addition to the edge shear, the plate is subjected to the compressive and tensile stresses in regions located opposite to and behind the shear loads, respectively. For this reason, the deformation waves make angles greater than 45 degrees with the x-axis and the clamped supports have led to local post-buckling deformation patterns. As Fig. 22 shows, apart from the edges, more kinks are induced in the interior regions of the plate, due to presence of the clamped supports and the resulting local deformations. Comparing Fig. 20 with Fig. 7 reveals that in the present case, in contrast to the simply supported plate, the buckling load of the pure metallic plate is not higher than that of the FGM plate. Moreover, buckling strengths of the present plate are lower than those of the simply supported plate (approximately, less than half), due to the compressive stresses exerted by the clamped edges.

21

8 7 6

PN

5 4 3 2

Pure Metal n=1 n=3 n=5

1 0 0

0.2

0.4 WN

0.6

0.8

Fig. 20 Shear post-buckling paths of a square FGM plate with clamped edges at x=0 and a, based on observations of the (x= a/2, y=b/2, z=h/2) reference point.

8 7 6

PN

5 4 3 2

Pure Metal n=1 n=3 n=5

1 0 0

0.5

1

1.5

WN Fig. 21 Shear post-buckling paths of a square FGM plate with clamped edges at x=0 and a, based on observations of the (x=a/4, y=b, z=h/2) reference point. Effects of the aspect ratio on the same plate (whose edges are clamped at x=0 and a) are studied in Figs. 23 and 24 that are plotted for a plate with a/b=2 aspect ratio. Comparing Figs. 10 and 23 show that present buckling strengths are slightly higher than one tenth of those of a simply supported plate. The

22

transverse post-buckling deformations of mid-plane of the plate are shown in Fig. 24. Since angle of the deformation region is larger than that of the simply supported plate, the two local buckling regions have constituted a single deformation band.

1

1.4 1.2

0.8

1 0.6 y (m)

0.8 0.6

0.4

0.4 0.2

0 0

0.2 0 0.2

0.4 0.6 x (m)

0.8

1

Fig. 22 Post-buckling transverse deformations of the mid-plane of the square FGM plate with clamped edges at x=0 and a (n=1 and PN = 6.62 ).

4 3.5 3

PN

2.5 2 1.5 1

Pure Metal n=1 n=3 n=5

0.5 0 0

0.5

1

1.5

2

2.5

3

3.5

WN Fig. 23 Shear post-buckling paths of a rectangular (a/b=2) FGM plate with clamped edges at x=0 and a, based on observations of the (x= a/4, y=b, z=h/2) reference point.

23

1 3 0.8

2.5 2

y (m)

0.6

1.5 0.4

1 0.5

0.2

0 0 0

0.5

1 x (m)

1.5

2

Fig. 24 Post-buckling transverse deformations of the mid-plane of the rectangular (a/b=2) FGM plate with clamped edges at x=0 and a (n=1 and PN = 3.97 ). 5.4 Influence of the elastic foundation on the shear post-buckling stiffness and deformations Influence of the elastic foundation on the shear post-buckling deformations and strengths are investigated in the present section for both square and rectangular FGM plates with simply supported edges. Stiffness of the elastic foundation has to be in the order of the stiffness of the plate to show pronounced effects. For this reason, a plate with the material information mentioned in Sec. 5.1 and the following specifications is considered: a / b = 1or 2, h / a = 0.025 or 0.05, n = 1 In order to plot the load-deflection curves, three reference points are chosen: (x= a/2, y=b/2, z=h/2), (x= a/4, y=b/2, z=h/2), and (x= a/5, y=4b/5, z=h/2) and the relevant post-buckling curves are plotted in Figs. 25 to 27, respectively. In these figures, results correspond to the Kw=0, 5e9, 1e10, and 5e10 N/m3 are compared. Although all these figures predict identical buckling loads for identical stiffnesses of the elastic foundation, Fig. 26 shows a notable change in the post-buckling curve for Kw= 5e10 N/m3; a fact that emphasizes improperness of the second reference point for tracing all the events. To confirm this evidence, post-buckling deformation patterns of the mid-plane of the FGM plate are depicted in Figs. 28 and 29 for Kw=0 and 5e10 N/m3, respectively, for an identical dimensionless shear load. Figs. 27 and 28 show that while the second reference point is located on the main post-buckling deformation band for the plate without elastic foundation, it is not so for a plate on an elastic foundation with Kw= 5e10 N/m3stiffness. Through comparing Figs. 25 to 27, it may be readily deduced that observations of the third reference point are the most consistent ones as this point remains on the deformation band. Comparing Figs. 28 and 29 lead to the conclusion that presence of the elastic foundation prevents formation of the global buckling bands; so that the resulting post-buckling deformation regions are local and their number is larger.

24

15

PN

10

Kw=0

5

Kw=5e9 Kw=10e9 Kw=50e9 0 0

0.5

1 WN

1.5

2

Fig. 25 Influence of stiffness of the elastic foundation on the post-buckling paths of the square FGM plate (n=1), based on observations of the (x= a/2, y= b/2, z=h/2) reference point.

15

PN

10

5

Kw=0 Kw=5e9 Kw=10e9

0

Kw=50e9

0

0.2

0.4

0.6

WN Fig. 26 Influence of stiffness of the elastic foundation on the post-buckling paths of the square FGM plate (n=1), based on observations of the (x= a/4, y=b/2, z=h/2) reference point.

25

15

PN

10

5

Kw=0 Kw=5e9 Kw=10e9

0 0

Kw=50e9

0.5

1

1.5

2

2.5

WN Fig. 27 Influence of stiffness of the elastic foundation on the post-buckling paths of the square FGM plate (n=1), based on observations of the (x= a/5, y=4b/5, z=h/2) reference point.

1 2

0.8

1.5

y (m)

0.6

1

0.4

0.5

0.2

0 0

0 0.2

0.4 0.6 x (m)

0.8

1

Fig. 28 Post-buckling transverse deformations of the mid-plane of the square FGM plate without elastic foundation (n=1 and PN = 13.91 ).

26

1 1.5

0.8

1

0.4

0.5

y (m)

0.6

0.2 0 0 0

0.2

0.4 0.6 x (m)

0.8

1

Fig. 29 Post-buckling transverse deformations of the mid-plane of the square FGM plate with elastic foundation (n=1 and K w = 5e10 N / m 3 , PN = 13.91 ). Post-buckling paths of the rectangular (a/b=2) FGM plate are traced using two distinct reference points: (x= a/2, y=b/2, z=h/2) and (x= a/4, y=b/2, z=h/2). The resulting post-buckling curves are plotted in Figs. 30 and 31, respectively, for various stiffnesses of the elastic foundation. Although both figures predict the same buckling loads for the identical stiffnesses of the elastic foundation, the post-buckling curves are slightly different. The post-buckling deformation pattern shown in Fig. 321, for a specific nondimensional shear load and Kw= 1e10 N/m3reveals that results of the first reference point may be more accurate. Moreover, comparing Fig. 32 with Figs. 18 and 19 leads to the conclusion that elastic foundation has increased number of the deformation bands from 3 to 5. Comparing results presented in Fig. 30 with those of Fig. 27 reveals that post-buckling strength of the rectangular plate is significantly lower than that of the square plate (when the b dimension is retained).

27

7 6 5

PN

4 3 2 1

Kw=0 Kw=1e9 Kw=5e9 Kw=10e9

0 -1.5

-1

-0.5

0

WN Fig. 30 Influence of stiffness of the elastic foundation on the post-buckling paths of the rectangular (a/b=2) FGM plate (n=1), based on observations of the (x=a/2, y= b/2, z=h/2) reference point.

7 6 5 PN

4 3 Kw=0

2

Kw=1e9 Kw=5e9

1

Kw=10e9

0 0

0.2

0.4

0.6 WN

0.8

1

1.2

Fig. 31 Influence of stiffness of the elastic foundation on the post-buckling paths of the rectangular (a/b=2) FGM plate (n=1), based on observations of the (x=a/4, y= b/2, z=h/2) reference point.

28

1 0.4

y (m)

0.8

0.3

0.6

0.2 0.1

0.4

0 0.2 0 0

-0.1 -0.2 0.5

1 x (m)

1.5

2

Fig. 32 Post-buckling transverse deformations of the mid-plane of the rectangular FGM plate ( n=1 and K w = 1e10 N / m 3 , PN = 5.63 ). 6 Conclusions In the present paper, shear post-buckling behavior of the FGM rectangular plates with movable or clamped edges is studied, using the most general 3D theory of elasticity and the full Green’s straindisplacement tensor. The finite element form of the governing equations is obtained by using a full consistent 3D Hermitian element with up to the second-order mixed derivatives continuity and solved by means of an arc-length solution method in conjunction with the modified Newton-Raphson procedure. Furthermore, effects of the elastic foundation on the post-buckling deformation/distortion shapes and strengths are investigated as well. Some of the practical conclusions drawn based on the results section may be summarized as follows: • For a square FGM plate, the post-buckling deformation wave makes approximately a 45 degrees angle with respect to the edges. • Proper choice of the reference point for measurement of the lateral deflection of the load-deflection (post-buckling) curve may be chosen only after studying the distribution of the transverse deformation of the entire plate for various magnitudes of the shear load. • Number of the buckling waves and distortions may increase by increasing the load magnitude (mode jumping). • Angles of the deformation waves change by changing the aspect ratio, load magnitude, and especially, stiffness of the elastic foundation. • The elastic foundation leads to local deformations and larger number of the deformation regions. In this case, choice of the reference point has to be accomplished with a special care. • Buckling strengths reduce remarkably when some of the edges are clamped. References [1] Lindberg HE, Florence AL. Dynamic pulse buckling: Theory and experiment. Dordrecht, The Netherlands: Martinus Nijhoff Pub.; 1987. [2] Donnell LH. Beam, plates, and shells. US: McGraw-Hill Inc.; 1976.

29

[3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26]

Jones RM. Buckling of bars, plates, and shells. Virginia: Bull Ridge Publishing; 2006. Budiansky B, Connor RW. Buckling stresses of clamped rectangular flat plates in shear. NACA Tech Note No. 1559; 1948. Smith ST, Bradford MA, Oehlers DJ. Elastic buckling of unilaterally constrained rectangular plates in pure shear. Eng Struct 1999;21:443–453. Azikov NS, Vasiliev VV, Paterekas AD. Buckling of composite plates under compression and shear. Mech Compos Mater 1990;2;351-3. Loughlan J. The shear buckling behaviour of thin composite plates with particular reference to the effects of bend-twist coupling. Int J Mech Sci 2001;43:771-792. Lopatin AV, Korbut YB. Buckling of clamped orthotropic plate in shear. Compos Struct 2006;76:94-98. Papadopoulos L, Kassapoglou C. Shear buckling of rectangular composite plates composed of concentric layups. Compos Part A 2007;38:1425–1430. Shufrin I, Eisenberger M. Shear buckling of thin plates with constant in-plane stresses. Int J Struct Stabil Dyn 2007;7:179–192. Shufrin I, Rabinovitch O, Eisenberger M. Buckling of laminated plates with general boundary conditions under combined compression, tension, and shear—A semi-analytical solution. ThinWalled Struct 2008;46:925–938. Wu L-Y, Wu C-H, Huang H-H. Shear buckling of thin plates using the spline collocation method. Int J Struct Stabil Dyn 2008; 8: 645–664. Tarján G, Sapkás Á, Kollár LP. Stability analysis of long composite plates with restrained edges subjected to shear and linearly varying loads. J Reinf Plast Compos 2010;29:1386-1398. Chirica I, Beznea E-F. Buckling behavior of the multiple delaminated composite plates under shear and axial compression. Comput Mater Sci 2012;64:173-178. Yang D, Huang Y, Li G. Shear buckling analysis of anisotropic rectangular plates. Yingyong Lixue Xuebao/Chinese J. Appl. Mech. 29 (2012) 220-224. Uymaz B, Aydogdu M. Three dimensional shear buckling of FG plates with various boundary conditions. Compos Struct 2013;96:670–682. Kosteletos S. Postbuckling response of laminated plates under shear load. Compos Struct 1992;20:137–145. Mittelstedt C, Erdmann H, Schröder KU. Postbuckling of imperfect rectangular composite plates under inplane shear closed-form approximate solutions. Arch Appl Mech 2011;81:1409-1426. Singh S, Kumar A, Postbuckling response and failure of symmetric laminates under in-plane shear. Compos. Sci Tech 1998;58:1949–1960. Singh S, Kumar A. Postbuckling response and strength of laminates under combined in-plane loads, Compos. Sci Tech 1999:59:727–736. Kim KD, Lomboy GR, Han S-C. A co-rotational 8-node assumed strain shell element for postbuckling analysis of laminated composite plates and shells. Comput Mech 2003;30:330–342. Han S-C, Lee S-Y, Rus G. Postbuckling analysis of laminated composite plates subjected to the combination of in-plane shear, compression and lateral loading. Int J Solids Struct 2006;43:57135735. Gupta A, Patel B., Nath Y., Postbuckling response of composite laminated plates with evolving damage. Int J Damage Mech 2014;23:222-244. Chen Q, Qiao P. Post-buckling behavior of imperfect laminated composite plates with rotationallyrestrained edges. Compos Struct 2015;125:117–126. Nie K, Liu Y, Dai Y. Closed-form solution for the postbuckling behavior of long unsymmetrical rotationally-restrained laminated composite plates under inplane shear. Compos Struct 2015;122:31–40. Raju G, Wua Z, Weaver PM. Buckling and postbuckling of variable angle tow composite plates under in-plane shear loading. Int J Solids Struct 2015;58:270–287.

30

[27] Ferreira AJM, Roque CMC, Neves AMA, Jorge RMN, Soares CMM. Analysis of plates on Pasternak foundations by radial basis functions. Comput Mech 2010;46:791-803. [28] Tornabene F, Fantuzzi N, Viola E, Reddy JN. Winkler–Pasternak foundation effect on the static and dynamic analyses of laminated doubly-curved and degenerate shells and panels. Compos Part B 2014;57:269–296. [29] Tornabene F, Fantuzzi N, Viola E. Interlaminar stress recovery procedure for doubly-curved, singly-curved, revolution shells with variable radii of curvature and plates using generalized higher-order theories and the local GDQ method. Mech Adv Mater Struc 2016;23:1019-1045. [30] Neves AMA, Ferreira AJM, Carrera E, Cinefra M, Jorge RMN, Soares CMM. Buckling analysis of sandwich plates with functionally graded skins using a new quasi-3D hyperbolic sine shear deformation theory and collocation with radial basis functions. ZAMM 2012;92:749-766. [31] Sofiyev AH, Kuruoglu N. Torsional vibration and buckling of the cylindrical shell with functionally graded coatings surrounded by an elastic medium. Compos Part B 2013;45:1133-1142. [32] Asemi K, Salehi M, Akhlaghi M. Three dimensional graded finite element elasticity shear buckling analysis of FGM annular sector plates. Aerospace Sci Tech 2015;43:1-13. [33] Shariyat M, Eslami MR. Dynamic buckling and postbuckling of imperfect orthotropic cylindrical shells under mechanical and thermal loads, based on the three dimensional theory of elasticity. Trans ASME, J Appl Mech 1999;66:476-484. [34] Alipour MM, Shariyat M. Semi-analytical buckling analysis of heterogeneous variable thickness viscoelastic circular plates on elastic foundations. Mech Res Commun 2011;38:594-601. [35] Mansouri MH, Shariyat M. Biaxial thermo-mechanical buckling of orthotropic auxetic FGM plates with temperature and moisture dependent material properties on elastic foundations. Composites: Part B 2015;83:88–104. [36] Shariyat M, Asemi K. Three-dimensional non-linear elasticity-based 3D cubic B-spline finite element shear buckling analysis of rectangular orthotropic FGM plates surrounded by elastic foundations. Compos Part B 2014;56:934-947. [37] Asemi K, Shariyat M, Salehi M, Ashrafi H. A full compatible three-dimensional elasticity element for buckling analysis of FGM rectangular plates subjected to various combinations of biaxial normal and shear loads. Finite Elem Analysis Des 2013;74:9–21. [38] Asemi K, Shariyat M. Highly accurate nonlinear three-dimensional finite element elasticity approach for biaxial buckling of rectangular anisotropic FGM plates with general orthotropy directions. Compos Struct 2013;106:235-249. [39] Shariyat M, Asemi K. 3D B-spline finite element nonlinear elasticity buckling analysis of rectangular FGM plates under non-uniform edge loads, using a micromechanical model. Compos Struct 2014;112:397–408. [40] Asemi K, Shariyat M. Three-dimensional biaxial post-buckling analysis of heterogeneous auxetic rectangular plates on elastic foundations by new criteria. Comput Meth Appl Mech Eng 2016;302:1-26. [41] Shariyat M, Asemi K. Uniaxial and biaxial post-buckling behaviors of longitudinally graded rectangular plates on elastic foundations according to the 3D theory of elasticity. Compos Struct 2016;142: 57–70. [42] Huang H-C. Static and Dynamic Analyses of Plates and Shells: Theory, Software and Applications. Springer-Verlag Berlin Heidelberg, 1989. [43] Reese S, Wriggers P, Reddy BD. A new locking-free brick element technique for large deformation problems in elasticity. Comput Struct 75 (2000) 291-304. [44] Yang J, Huang Q, Hu H, Giunta G, Belouettar S, Potier-Ferry M. A new family of finite elements for wrinkling analysis of thin films on compliant substrates. Compos Struct 2015;119:568–577. [45] Hazrati H, Shariyat M. Natural frequency & shape mode of isotropic and FGM plates using firstorder shear deformation theory with finite element method and investigation of rotation of the line of nodes of vibration. The 2nd International Conference on Composites: Characterization, Fabrication and Application (CCFA-2), Kish Island, Iran, December 27-30, 2010.

31

[46] Zienkiewicz OC, Taylor RL. The Finite Element Method, Volume 1: The Basis. ButterworthHeinemann, Oxford, 2000. [47] Eslami MR. Thermo-mechanical buckling of composite plates and shells. Tehran: Amirk-abir University Press; 2010. [48] Shen H-S. Functionally graded materials: Nonlinear analysis of plates and shells. CRC Press, Taylor & Francis Group; 2009. [49] Riks E. The application of Newton’s method to the problem of elastic stability. J Appl Mech 1972;39:1060-1066. [50] Wempner GA. Discrete approximations related to non-linear theories of solids. Int J Solids Struct 1971;7:1581-1599. [51] Crisfield MA. A fast incremental/iterative solution procedure that handles “snap through”’. Comp Struct 1981;13:55-62. [52] Ramm E. Strategies for tracing the nonlinear response near limit points. In: Non-linear Finite Element Analysis in Structural Mechanics (Eds. Wunderlich W, Stein E, Bathe KJ), Berlin: Springer-Verlag, 1981; 68-89. [53] Crisfield MA. Non-linear finite element analysis of solids and structures. Chichester, England: John Wiley & Sons Ltd.; 2000. [54] Carrera E. A study on arc-length-type methods and their operation failures illustrated by a simple model. Comput Struct 1994;50:217-229. [55] Hellweg H-B, Crisfield MA. A new arc-length method for handling sharp snap-backs. Comput Struct 1998;66:705-709. [56] Kweont JH, Hong CS. An improved arc-length method for post-buckling analysis of composite cylindrical panels. Comput Struct 1994;53:541-549. [57] Zhu J-F, Chu X-F. An improved arc-length method and application in the post-buckling analysis for composite structures. Appl Math Mech (English Edition) 2002;23:1081-1088.

32