Static deﬂection analysis of ﬂexural rectangular micro-plate using higher continuity ﬁnite-element method

– In this paper, strain gradient theory is used in developing a mathematical model based on classical ﬂexural Kirchhoﬀ plate theory that can predict static response of rectangular micro-plates. The result of this new model is a sixth order diﬀerential equation. Order of diﬀerential terms in Galerkin weak form of the equation is reduced so that C 2 hierarchical p -version ﬁnite elements with second order global smoothness can be used to solve the problem. With diﬀerent boundary conditions, the computed deﬂection distribution of micro-plates is compared with those of the classical theory, in which length scale parameters are not present. A series of studies have revealed that when length scale parameters are considered, deﬂection of a rectangular plate decreases with increasing the length scale eﬀect; in other words micro plates exhibit more rigidity than what is predicted by the classic model. Here, deﬂections are normalized with respect to results obtained from classical plate theory. Comparison of maximum deﬂection values obtained from the extended model for micro plates with those available from the classic plate model indicates that classical theory overestimates displacement values and the largest error is observed for square micro plates. The overestimation levels oﬀ for plates with aspect ratios greater than three.


Introduction
In classical theory of linear elasticity, strain is defined based on local distribution of the stress field.In this model nonlocal effects and length scale parameters are not considered, thus rendering the classical plate theory as inappropriate for describing the mechanical behavior of a micro-structure, e.g.polymeric foams, high-toughness ceramics, high strength metal alloys, granular materials or porous bones.Behavior of such structures is strongly influenced by nonlocal effects in the stress field and by significance of the internal length scale effects.Experimental studies have demonstrated the effects of microstructure on mechanical behavior of micro-plates [1][2][3][4][5][6][7].Microstructural effects become important when structures, such as bars, beams or plates have extremely small overall dimensions with respect to the internal length scale of their constitutive material.Structures of this size find applications in microelectromechanical systems (MEMS) or nanoelectromechanical systems (NEMS).The load response of these structures cannot be described properly using the classic linear elasticity model.For the mentioned cases of a Corresponding author: a.ahmadi@kgut.ac.ir materials and structures, one has to use generalized or higher-order theories of linear elasticity; which are characterized by non-locality of the stress field.A comprehensive physical description of the length scale parameter is still under discussion, nonetheless the concept has been successfully applied in modeling a variety of problems in microstructures.Results obtained from higher order models are also comparable to those obtained from various homogenization techniques [8].Mindlin [9] and Mindlin and Eshel [10] generalized the elastic theory of microstructured solids.Eringen [11] developed the theory of micropolar elastic solids in accordance to the Cosserats theory of elasticity [12]; and couple stress theory is due to Toupin [13] and Koiter [14].The general theory of Mindlin [9] includes three equivalent forms which are defined on the basis of three different expressions for strain energy density (the first expression involves gradients of displacements, the second, gradients of strain and the third, gradients of rotation) and also three length scale parameters.Couple stress theory is based on the third expression of strain energy density.Simplified strain gradient theory of Mindlin and Eshel [10] for the case of an isotropic linear elastic microstructured solid involves five elastic constants in addition to the two classical Lamé's constants and contains only first-order strain gradients in the strain energy density expression.A review of the above-mentioned higher-order theories of elasticity can be found in the works of Tiersten and Bleustein [15], Vardoulakis and Sulem [16], Lakes [17].Papargyri-Beskou and Beskos [18] investigated the static flexural response of rectangular micro-plates.Also, Park and Gao [19] obtained a closed form solution of a simple shear problem.In this paper, static deflection of a rectangular microplate, including one length scale parameter, is investigated by utilizing C 2 hierarchical p-version finite-element setting.Effects of length scale parameter are then studied for SSSS, CCCC, CCSS, CCF F and SSF F boundary conditions of rectangular micro-plate and relevant observations are discussed.Moreover, micro-plates with different aspect ratios for several boundary conditions of microplate are studied in detail.

Strain gradient elasticity
In order to account for non-local effects that are present in micro-structures, the relation between stress and strain fields can be generalized as follows [17,20]: where σ 0 is the initial stress, ε is the strain variable, g is length scale parameter and η is strain gradient.Simplest possible version of the strain gradient elasticity theory due to Mindlin [9,10] (which has five constants besides two Lamé's constants) is a model with just one constant in addition to the Lamé's constants.The constitutive equations for this simple model are given as [18]: Here, σ and τ are the total and the classical Cauchy stress tensors respectively, Ĩ is the unit stress tensor; ε and tr(ε) are the strain tensor and its trace which are expressed in terms of the displacement vector u as: Parameter g 2 is the volumetric strain energy gradient coefficient or simply gradient coefficient with g being the internal or characteristic length of microstructure grains; λ and μ are the two classic Lamé's constants.
Imposing the strain gradient term with length scale parameter on the conventional elasticity model as a constraint was first discussed by Farahmand and Arabnejad [20].Comparison of experimental results obtained from torsion and bending tests of beams with theoretical ones obtained from this and other higher-order elasticity models have revealed that magnitude of the gradient coefficient g (internal length) is of the same order as diameter of the basic building block in a microstructure, e.g., the grain in metals or ceramics, the osteon in bones or the cell in foams [16,17,[21][22][23].

Governing equation for flexural bending of micro plates
Based on the work by Papargyri-Beskou and Beskos [18], invariant form of the mathematical model in terms of transverse displacement field, w(x) of a microplate subjected to static loading is [18]: where plate's rigidity D = Eh 3  12(1−ν 2 ) .Weak formulation of Equation (7) yields where, v represents all admissible weight functions, which for Galerkin formulation becomes v = δw.Repeated application of Green's identity (integration by parts) yields a symmetric weak form for Equation (8) as follows: where the boundary integral terms resulting from application of Green's identity are:  In this study, all numerical experiments are conducted on rectangular micro plates with dimensions of a and b and thickness of h.Finite element approximation space employed to obtain the results presented herein is constructed from a 100 element uniform mesh of two dimensional C 2 p-version elements evaluated at p-level P ξ = P η = 5.
In studying micro-plates, boundary conditions are categorized into two main groups of Classical Boundary Conditions (C.B.C) and Non-Classical Boundary Conditions (N.B.C).The (N.B.C) for a simply supported micro plate are expressed as follows [24] w xx = w yy = w xxxx = w yyyy = w xxyy = 0 at x = 0, a w yy = w xx = w yyyy = w xxxx = w yyxx = 0 at y = 0, b (11) and for clamp boundary conditions (N.B.C) are obtained as follows [24] w xx = 0 or Moreover, the Boundary Integral Terms (B.I.T) also referred to natural boundary conditions can be divided into two groups of classical and non-classical terms.Nevertheless, the non-classical terms of (B.I.T) vanish along intersegmental boundaries through assembly and also disappear on exterior boundaries of the micro-plate, where essential boundary conditions are specified.Here, the non classical terms (N.B.I.T) are:

Finite element approximation
The weak form given by ( 8) requires second order global smoothness as a minimum for integration to be valid in the Lebesgue sense; i.e. the trial solution w ∈ W 3 , where the trial space W 3 is defined as: where A basis can be formed [25] for trial space W 3 through tensor product of spaces U 3 and V 3 where and A p-version element with 2nd order global smoothness was constructed in this manner to approximate the weak solution.
Starting with one dimensional p-version C 0 finiteelement (Fig. 1) which interpolates w as w = where p is the degree of interpolating polynomial and basis N k (ξ) are denoted as where - for i = 2, 3, . . .p.A C 2 p-version element as shown in (Fig. 2) can be generated by applying nodal conditions = φ 3i for i = 1, 2 so that the basis vectors are transformed from C 0 to C 2 .This transformation results in elimination of four basis functions from the C 0 element's hierarchical node.Degrees of freedom for C 2 element are ∂ i w(ξ) ∂ξ i where i = 0, 1, 2 at element end nodes and ∂ j w(ξ) ∂ξ j ξ=0 where j = 6, 7, . .., p at element's hierarchical node.Now, a two dimensional C 2 element can be generated through tensor product of two orthogonal one dimensional C 2 elements; with degree of polynomial p = P ξ for interpolation in ξ direction and p = P η in η direction, Figure 3.
The degrees of freedom at each corner point of this element are: (20) Degrees of freedom at hierarchical nodes are of the form ∂ k w ∂ξ i ∂η j .Basis functions generated for this element yield interpolants that are globally smooth up to second order derivatives; however, because of the tensor product operation, some terms from the third and fourth derivative sets appear at element corner nodes.Even though this element satisfies requirements of the trial space, W 3 , nonetheless, the actual geometry of this element in the physical domain is limited to rectangular shapes.The reason for that is the incomplete appearance of higher derivatives (i.e.third order and fourth) as degrees of freedom in corner points, which makes the transformation of nodal degrees of freedom from local ξ−η to global x−y space impossible.Nevertheless, this element is sufficiently smooth for investigating the effects of length scale on bending behavior of microstructure plates.
Entries of the element stiffness matrix [K] in terms of element basis functions are as follows: and entries of elements' load vector are given as are:

Numerical experiments
The C 2 elements used here are generated through tensor multiplication of two one-dimensional C 2 elements which are not difficult to generate.It is true that higher continuity elements do not exist in commercial softwares and their generation/implementation is not a straight forward matter, nonetheless they are needed to solve the problem at hand.The alternative to using a C 2 element would be to utilize C 0 elements with an equivalent differential formulation.In an equivalent formulation the sixth order differential equation has to be written as an equivalent system of three second order differential equations.After applying the weak form to such a system, theoretically one would be able to use C 0 elements to solve the system.This approach would triple the problem size immediately; now one has to deal with three variables instead of one.Such an approach is used in commercial packages to solve the fourth order differential equations for plates and shells; which give rise to numerical difficulties such as shear locking.In a similar manner, one should expect unseen numerical difficulties when dealing with the equivalent second order differential system.Based on our experience, solving the stronger system, i.e. sixth order differential equation using C 2 element would have less computational issues than solving three second order differential equations which could be cast in different forms [26][27][28].In order to verify the implemented finite element formulation and methodology, the results are compared to a closed form solution presented in reference [18].The solution has been reported for simply supported boundary conditions [18], for this solution the unknown displacement field and loading function are represented using Fourier basis.This type of consideration satisfies the boundary conditions [29] and the solution is found as follows [18] w Then substituting Equation ( 23) into Equation ( 6) yields The classical solution is found by letting g = 0 and displacement is obtained as: Denoting W max as the maximum value of W mn (Eq.( 24)) for n = 1, the value of m for which W mn becomes maximum is found by setting dW mn /dm = 0. Thus, for a square plate with side length a, and n = m = 1 one obtains Wmax W c max in the form (Ref. [18]) as, Figure 4 demonstrates the accuracy of the implemented finite element procedure here for studying the static response of micro plates.Figure 5 for the case with CCCC boundary conditions, drastic decrease in the nondimensionalized maximum deflection is observed; which is indicative of high stiffness due to the inclusion of length scale parameter.Moreover, minor changes in nondimensionalized maximum deflection are obtained for SSF F boundary condition.Study of micro-plates with free boundary conditions along the boundary with dimension b indicates that the nondimensionalized maximum deflection Wmax W c max rises steeply with increasing b/a ratio.All free cases, as to be anticipated, converge to micro beam solution.Also note that the rising rate increases with increasing length scale values.Noting from Figures 6-10 all nondimensionalized maximum deflection curves rise quickly for all cases and each reaches a steady value for dimension ratios b/a > 3.This feature reflects the diminishing effects of the constant-length boundaries for b/a > 3. Comparing all boundary conditions, CCCC boundary condition for micro plate shows the least amount of change due to variations in micro-plate's aspect ratio.

Conclusion
When internal length scale parameter of a microstructure is comparable to the length scale of its overall      geometry, classical elasticity model fails.Hence, in describing its mechanical response, higher-order theories are required to construct proper models.In this paper the simplest gradient elastic material model, i.e. one constant (internal length scale parameter) in addition to two Lamé's constants of the classical elasticity has been considered and studied.The resulting elastic model which is based on gradient elasticity and describes the steady state response of flexural Kirchhoff plates subjected to mechanical loading is a sixth order differential model instead of the usual fourth order plate bending model.Higher continuity p-version finite elements were utilized in solving the resulting sixth order differential equation in the weak sense.Based on numerical experiments, it is concluded that increasing the length scale parameter leads to decrease in deflection; the behavior is due to higher rigidity caused by the length scale effects.Variations in aspect ratio of the micro-plate have significant effects on rigidity, for b/a > 3, however, this effect is related to mere geometrical characteristic of plate and it is significant effect in all cases except in the case of an all clamped micro plate.

Fig. 4 .
Fig. 4. Comparing normalized deflection Wmax W c max variation for different values of length scale ratio ( g a ) 2obtained from finiteelement method with exact solution given in reference[18].

Fig. 5 .Fig. 6 .
Fig. 5. Normalized deflection Wmax W c max variation of a square micro plate versus length scale ratio ( g a ) under different boundary conditions.

Fig. 7 .
Fig. 7. Normalized deflection Wmax W c max variation of rectangular micro plates with different aspect ratios and for different values of length scale ratio ( g a ) under SSF F boundary condition.

Fig. 8 .
Fig. 8. Normalized deflection Wmax W c max variation of rectangular micro plates with different aspect ratios and for different values of length scale ratio ( g a ) under CCCC boundary condition.

Fig. 9 .
Fig. 9. Normalized deflection Wmax W c max variation of rectangular micro plates with different aspect ratio and for different values of length scale ratios ( g a ) under CCSS boundary condition.

Fig. 10 .
Fig. 10.Normalized deflection Wmax W c max variation of rectangular micro plates with different aspect ratios and for different values of length scale ratio ( g a ) under SSSS boundary condition.