## Services on Demand

## Journal

## Article

## Indicators

- Cited by SciELO
- Access statistics

## Related links

- Cited by Google
- Similars in SciELO
- Similars in Google

## Share

## Latin American Journal of Solids and Structures

##
*On-line version* ISSN 1679-7825

### Lat. Am. j. solids struct. vol.11 no.6 Rio de Janeiro Nov. 2014

#### https://doi.org/10.1590/S1679-78252014000600007

**An efficient formulation for linear and geometric non-linear membrane elements**

**Mohammad Rezaiee-Pajand ^{*}; Majid Yaghoobi**

Civil Engineering Department, Fer-dowsi University of Mashhad, Iran

**ABSTRACT**

Utilizing the straingradient notation process and the free formulation, an efficient way of constructing membrane elements will be proposed. This strategy can be utilized for linear and geometric non-linear problems. In the suggested formulation, the optimization constraints of insensitivity to distortion, rotational invariance and not having parasitic shear error are employed. In addition, the equilibrium equations will be established based on some constraints among the strain states. The authors' technique can easily separate the rigid body motions, and those belong to deformational motions. In this article, a novel triangular element, named SST10, is formulated. This element will be used in several plane problems having irregular mesh and complicated geometry with linear and geometrically nonlinear behavior. The numerical outcomes clearly demonstrate the efficiency of the new formulation.

**Keywords:** Strain states, Geometric non-linear formulation, Plane problems, Finite element.

**1 INTRODUCTION**

Various types of formulation have been established since 1980s. The two most famous ones are named free formulation (Bergan and Nygard, 1984) and assumed strain scheme (MacNeal, 1978). Based on these works, the high-performance elements were formulated. Some of the conspicuous characteristics of these elements are: simplicity, rapid convergence, yielding similar accuracy in displacements and strains, rotational invariance, rank sufficiency, insensitivity to geometric distortion and being mixable with other elements (Felippa, 1994).

During the early 90s, the parameterized variational principle dramatically changed the knowledge towards highly efficient element formulations. Based on this principle, researchers managed to define a continuous space of the elastic functional. In the early 1990s, the parameterized variational principle was applied in the finite element method. Consequently, the formulation of high-performance elements was developed more fully, and a continuous space of the elastic functionals was provided. Making the continuous space of the functionals stationary, produces free parameters for the formulation of the element (Felippa and Militello, 1990). As a result, this will lead to the construction of finite element templates. Note that the optimization process of templates is quite complex and requires innovations. The large number of the free parameters, symbolic processing and optimizing the entries of the matrix are the difficulties which are encountered in the process.

As Bergan and Nygard were attempting to develop free formulation, John Dow et al. proposed another scheme named strain gradient notation (Dow, 1999). In the latter technique, the displacement interpolation field is expressed in terms of strain states. In this way, the roots of many errors in finite element method can be investigated and eliminated. It should be added that strain gradient notation is an explicit and simple representation of free formulation. Utilizing the strain gradient notation process and the free formulation, displacement and strain fields were written in terms of strain states. In addition, several optimality conditions were proposed. By using these optimality criteria, the resultant element becomes rotational invariance and insensitive to geometry distortion. Moreover, the parasitic shear error is eliminated, and the equilibrium equations are satisfied in this formulation.

Co-rotational technique is usually employed in finite element formulation for the geometrically nonlinear analysis of structures. Although, this scheme has been widely utilized in construction of solid (Norachan et al, 2012), shell (Li and Vu-Quoc, 2007; Battini, 2007), bending plate (Izzuddin, 2005) and beam elements (Battini, 2008a; Urthaler and Reddy, 2005), it has been rarely employed in plane structures. In this tactic, the rigid body motions are separated from deformations causing strains. This goal is achieved by attaching a local reference frame to the element. This coordinate system rotates and translates with the elements of the structure. To relate the global coordinate system to the co-rotating one, a transformation matrix is required. It should be added that some of the formulations in linear space are by themselves complicated. Hence, it is quite difficult to transform them from a linear space to a nonlinear one. This is particularly the case in which the template process proposed by Felippa is used (Felippa, 2006; Battini, 2008b). There is no need for transforming matrix in the mentioned formulation, and a general coordinate can realize the rigid body motions.

In the following sections, the simplicity and robustness of suggested technique in the geometrical nonlinear analysis of structures will be clearly demonstrated. Based on proposed formulation, a new triangular element, named SST10, is suggested. It is worth emphasizing that the suggested element is useful to analyze the plane problems with irregular mesh and complicated geometry. Note that the triangular elements are able to model the geometry of the structure more properly. However, the low-order triangular elements are less accurate than the quadrilateral elements (Eom et al., 2009). It should be added that various methods have been deployed to construct triangular elements (Felippa, 2003a; Zhou and Sze, 2012; Pan and Wheel, 2011; Huang et al., 2010; Caylak and Mahnken, 2011). In the following sections; several numerical tests are performed to demonstrate the capabilities of the proposed formulation. The results indicate high accuracy, low sensitivity to geometry distortion, rotational invariance and rapid convergence. Furthermore, it is proven that SST10 element performs properly in nonlinear numerical tests.

**2 STRAIN GRADIENT NOTATION FORMULATION**

In strain gradient notation, Taylor series expansion of the strain field is expressed around the coordinates' origin (Dow, 1999). In the two-dimensional state, the strain field has three coming appearances:

In these equations, o as a subscript denotes the gradient of strains in the origin. It should be reminded that magnitudes of the strain gradient at the origin are called the strain states. (e_{x})_{o} is the magnitude of axial strain *e _{x} * at the origin,

*(e*)

_{x,x}_{o}and

*(e*

_{x,y})

_{o}are the rate of

*e*variation in

_{x}*x*and

*y*directions, in the origin's vicinity, respectively. Similarly, other coefficients can be defined. To obtain the displacement field, the coming strain-displacement equations and rotation function are applied.

By considering the strains and rotations in the origin, linear parts of the displacement field can be determined. It is important to note that coefficients of the higher-order terms of the displacement field can be similarly obtained by deploying derivatives of the strains. Finally, the element's fields can be expressed as below:

Where M *u*_{o}and *v*_{o} are the rigid translations in *x * and *y * directions, respectively.

**3 OPTIMALITY CONDITIONS**

To construct the efficient elements, the optimization condition can be inserted into the formulation. In this way, the roots of many errors may be recognized. In the following sections, several optimali-ty constraints are introduced.

**3.1 Satisfying the equilibrium equations**

For a homogeneous elastic continuum in a plane-stress or plane-strain state, the equilibrium equations are written in the following form:

Where, the functions *F _{x }(x,y)* and

*F*are respectively the force field functions along the

_{y}(x,y)*x*and

*y*directions in a rectangular coordinate system. It is worth reminding that in planar problems, the gradient of the force field, in the perpendicular direction

*(z)*to the plane of the element, is equal to zero. The stress fields are defined in the rectangular coordinate system. Utilizing the stress-strain relations, equation 5 can be rewritten in terms of the strain fields. Considering Hooke's law for a homogeneous elastic condition, the following relations are valid for a plane problem:

In the above equation, *l* will be equal to , and for the plane-strain and plane-stress cases, respectively. The shear modulus, Poisson's ratio and elastic modulus, are denoted by *G, v *and *E, *correspondingly. Introducing equation 6 into relation 5 will result in the following equilibrium equations:

**3.2 Planar pure bending test**

This test has been utilized for finding the optimum bending templates by Flippa (Felippa, 2003a, 2003b, 2006). In this work, Euler-Bernoulli beam has been applied to assess the response of the template to planar bending in *x * and *y * directions. For this purpose, the energy ratios (r) are measured. It should be added that *r _{x} * and

*r*denote the bending energy ratios of the beam's rectangular part in

_{y}*x*and

*y*directions, respectively. Note that an element will be capable of modeling the bending in any arbitrary directions when r

_{x}= 1 or

*r*= 1 (r = 1). On the other hand, an element will be over stiff or over flexible while

_{y}*r >*1 or

*r <*1. Any aspect ratio of the element will have optimum behavior in bending when (r = 1). If (

*r*) increases as the aspect ratio soars, the locking shear will appear in the element.

At this stage, the mentioned test is investigated based on the strain states. In planar bending in *x * direction, the real stress field varies linearly along the *y * direction. By employing the Hook's law, the strain field can be given as the next form:

Where, *b _{1}, b_{2}, b_{3} * and

*b*are constant coefficients. According to equations 8, it is obvious that only(

_{4}*e*)

_{x}_{o},(

*e*)

_{y}_{0}

*,(e*and (

_{x})_{o}*e*)

_{y,y}_{o}participate in the aforementioned strain field. Similarly, for bending in

*y*direction,

*(e*)

_{x}_{o}, (

*e*)

_{y}_{0},

*(e*)

_{x}_{x}_{0}and (

*e*)

_{y,x}_{o}are incorporated into the real strain field. Hence,

*(e*, (

_{x})_{o}*e*)

_{y}_{0}, (

*e*)

_{x,x}_{0}, (

*e*)

_{x,y}_{0}, (

*e*)

_{y,x}_{0}and (

*e*) required for achieving the real responses in planar bending case. In addition, the element assumed strain field should always include constant strain states and the rigid motions states. This condition is necessary for convergence. It is worth emphasizing that the planar bending test based on the strain states can be used for any element with any geometries and meshes. In other words, not only it can be utilized for triangular and rectangular elements, but also it can be employed for other elements.

_{y,y}**3.3 Rotational invariance**

In a number of the specific elements, properties of the element will not change with rotation of the coordinate system. These elements fall in the category of the rotational invariance elements. Since elements can be configured with different rotational orientation in the mesh of a structure, the element should be rotational invariance. In order to have the latter property, terms of the strain field with a complete order should be considered (Dow, 1999). For example, rotational mapping of a constant strain state is as follows:

According to this equation, an element is capable of representing constant strains with respect to any system of the coordinates, only if its formulation takes into account all three cases of the constant strain states. It should be mentioned that invariance constraint has been used by researchers in the development of fine element templates (Felippa, 2000; Felippa and Militello, 1999).

**3.4 Lack of the parasitic shear error**

The appearance of axial strain states in the shear strain interpolation function causes the parasitic shear effect, which itself will lead to stiffen of the element (Dow, 1999). In order to gain knowledge about parasitic shear error, the formulation of the four-node rectangular element is examined. The interpolation functions of this element consist of an incomplete polynomial which contains neither *x ^{2 }* nor

*y*terms. The strain gradient notation of these terms is as follows:

^{2}Based on the strain-deformation relations (equations 2), the interpolation functions of the element can be written in the following form:

The advantages and the deficiencies of this model can be investigated through Taylor's expansion of this polynomial. By examining the shear strain series, it is noticed that the shear strains are independent of the axial strain. In the shear strain interpolation function, two strain states (*e _{x,y}*)

_{o}and (

*e*)

_{y,x}_{o}improperly appear. Hence, if the element undergoes flexural deformations, the strain states will be nonzero and will incorrectly represent a portion of the shear strain.

For elements which are formulated using the strain gradient notation, the parasitic shear error can easily be eliminated by excluding the incorrect strain states from the shear strain polynomial. The parasitic shear error decreases as the finer meshes are used. However, even coarser meshes can produce good results if elements are exempted from this error. It should be noted that utilizing complete interpolation functions will prevent the appearance of this error. Considering the abovemen-tioned statements, the following strain states will be used for the model:

**4 LINEAR FORMULATION**

Using the abovementioned strain states (equation 12), the displacement field of this model is defined as below:

The strain field is obtained as follows:

It should be added that *x-y * coordinates is a general one, and can be assumed anywhere in the problem domain. Based on relation 14, equation 7 can be written as below:

In these relationships, the body forces are assumed to be zero. According to equations 15, both strain states, (g_{Xy},_{x})_{o} and *(g _{xy,y})_{0} *can be expressed in terms of other strain states. By employing these equations, the number of the unknowns will be decreased to 10. Furthermore, the vector of the strain states can be given in the subsequent shape:

According to this equation, the strain interpolation field is written in the below matrix form:

In addition, the displacement interpolation field is obtained in matrix form as follows:

Inserting nodal coordinates into relation 19 yields the coming result:

In which **D** denotes the displacement vector. The matrix is obtained by placing the nodal coordinates into . It should be mentioned that provides a relation between the strain state vector and the displacement vector. Therefore, the boundary conditions can be applied to the formulation. In addition, the strain and displacement interpolation fields can be expressed in the succeeding form:

On the other hand, the potential energy function is attained in terms of the elastic matrix (**E**) as follows:

By utilizing equation 25, the stiffness matrix can be easily computed. The stiffness matrix (**K**) is also obtained as follows:

**5 GEOMETRICALLY NONLINEAR FORMULATION**

In this study, the geometrically nonlinear formulation is carried out by employing a co-rotational tactic. In this way, the displacement function corresponding to the rigid motions can be written in the succeeding shape:

In these relationships, *u*_{o}and *v*_{o} denote the rigid translations in *x * and *y * directions, respectively. Moreover, the rigid rotation of the coordinates' origin is shown by *q _{o}. * By transferring the origin of the x-y general coordinates to center of the element, the local coordinates will be in hand. This coordinates will be established in the initial shape and will be unchanged during the deformations. Having

*u*'

_{o}*v*

_{o}and

*q*

_{o}for each element, the rigid body motions will be found by equations 27. By utilizing relation 13, the deformations causing strains can be calculated with help of the next equalities:

Consequently, displacements of the element can be expressed in the local coordinates system as follows:

To find the strains, equations 28 can be used. In fact, the strain field similar to relations 14 will be obtained. In linear formulation with small strain, *q*_{o} is tiny, and it can be negligible. As a result, sin *q _{o} * and cos

*q*are assumed to be equal to (r

_{o}_{r})

_{o}and 1, respectively. In this case, the equalities 29 and 13 will be the same.

According to equations 15, both strain states, *(g _{xy,x})_{o }and (g_{xy,y})_{0} *can be expressed in terms of other strain states. By employing these equations, the number of the unknowns will be decreased to 10. Furthermore, the vector of the strain states can be given in the subsequent shape:

The strain and displacement interpolation field can be given in the succeeding matrix form:

By inserting the nodal coordinated into equation 33, the vector of nodal displacement can be obtained. According to equations 31 and 33, the interpolation field of the displacements and strains are written in terms of nodal displacement vector as follows:

Assuming*K*_{q} = *dv* yelds the coming result:

Where, *U *denotes the strain energy. By calculating the derivative of internal strain energy function with respect to nodal displacements, the internal load vector can be computed as the next form:

In equation 43, **P**_{1} is a vector which includes the first nine entries of the internal load vector and its tenth entry is incorporated into **P**_{2}. At this stage, by taking the derivative of internal load vector, the tangent stiffness matrix can be written in the coming shape:

In the last relations, **K**_{11} , **K**_{12} and **K**_{22}are 9×9, 9×1 and 1×1 matrices, respectively.

**6 THE ELEMENT'S DEGREES OF FREEDOM**

As previously mentioned, ten nodal unknowns should be identified to formulate the element. The geometry of the proposed element is shown in Figure 1.

In order to calculate* *and* *, the displacement in the direction perpendicular to the side of the element is expressed in terms of *u *and *v. *In Figure 2, this is carried out for one of the element's side shown by ij.

**7 NUMERICAL TESTS**

In the first stage, the proposed element is employed in several linear numerical analyses. The obtained results are compared with the responses of other researchers' well known elements. In the most cases, the triangular elements have not been used by other researchers. Therefore, for the sake of comparison, the responses of famous quadrilateral elements are deployed. Finally, the proposed element is applied to three geometrically nonlinear tests. In the following paragraphs, the characteristics of the elements applied in the numerical analysis are introduced:

1. Triangular element with rotational degrees of freedom; FF (αe=1.5, [β=0.5) (Bergan and Felippa, 1985).

2. Cook's plane hybrid triangle; CSTHybrid (Cook, 1987; Eom et al, 2009).

3. The modified enhanced assumed strain triangle; ME AS (Yeo and Lee, 1996; Eom et al., 2009).

4. Optimally fabricated assumed natural deviatoric strain triangle; OPT (Felippa, 2003a; Eom et al., 2009).

5. Airman's triangular element with spurious mode control; Allman (Allman, 1984; Eom et al., 2009).

6. Allman 88 triangular element integrated by 3-midpoint rule; Allman(3m) (Felippa, 2003a; All-man, 1984).

7. The triangular Hybrid Trefftz element with rotational degrees of freedom; HTD (Choo et al., 2006).

8. Retrofitted LST with ^{α}^{b=4/3}; LST-Ret (Felippa, 2003a).

9. Element with internal parameters and formulated by the QACM-I; AGQ6-II (Chen et al., 2004; Cen et al., 2009).

10. Stress hybrid element; PS (Pian and Sumihara, 1984; Chen et al., 2004; Cen et al., 2009).

11. Membrane element with drilling degrees of freedom; Q4S (Macneal and Harder, 1988; Cen et al., 2009).

12. Non-conforming isoparametric element with internal parameters; QM6 (Taylor et al., 1976; Chen et al., 2004; Choi et al., 2006; Cen et al, 2009).

13. Non-conforming isoparametric element with internal parameters; Q6 (Wilson et al., 1973; Cen et al., 2007, 2009).

14. Ibrahimbegovic's plane element with true rotation; IB (Ibrahimgovic et al., 1990; Choi et al, 2006).

15. The quadrilateral element with two enhanced strain modes; QE2 (Piltner and Taylor, 1995, 1997; Cen et al, 2009).

16. Assumed strain element; B-Q4E (Piltner and Taylor, 1997; Cen et al., 2009).

17. 4-node membrane elements with analytical element stiffness matrix; QAC-ATF4 (Cen et al., 2009).

18. A 4-node quadrilateral element with one-point quadrature integration procedure; Qnew (Fred-riksson and Ottosen, 2004).

19. Mixed four-node elements based on the Hu-Washizu functional; HW14-S and HW18-SS (Wisniewski and Turska, 2009)

**7.1 Linear tests**

In order to gain insight into the effectiveness of the proposed triangular element, seven linear tests, which are introduced by other researchers, are solved. These problems have formerly been used by the other researchers for testing a variety of different elements, and their results are available.

**7.1.1 Cantilever beam with five-element irregular mesh**

The geometry, loading and meshing of this structure are shown in Figure 3. The elasticity modulus and Poisson's ratio are 1500 and 0.25, respectively. All units of the problem are consistence. To analyze the beam, two load cases are applied: 1- bending moment *M, * which induces pure bending 2- concentrated force P, which produces linear bending. Under these load cases, the accurate vertical displacements of point *A * are 100.00 and 102.6, respectively (Cen et al., 2009). As previously mentioned, SST10 element is deployed to analyze the aforementioned beam. Figures 3.2 and 3.3 show two types of meshing for this structure, and the responses of these cases are calculated. Since other researchers have not used triangular elements for the analysis of this beam, the responses of their well-know quadrilateral elements are given, for the sake of comparison. The responses of different elements are listed in Table 1. It is obvious that SST10 element reaches to the exact solutions when the beam is subjected to the first load case. Under the second load case, the magnitude of error is approximately *2%. *It can be proven that the proposed element is insensitive to geometry distortion because the responses of meshes *A * and *B * are roughly the same.

**7.1.2 Cantilever beam with four-element irregular mesh**

Four irregular quadrilateral elements are utilized for mesh of the cantilever beam. In addition to this, the aforesaid beam is meshed in two different ways by using SST10 element, as shown in Figure 4. The elasticity modulus, Poisson's ratio and the beam's thickness are 30000, 0.25 and 1, respectively. Once again, all units of the problem are consistence. The geometry of the structure is shown in Figure 4.

At the end of the beam, a distributed parabolic shear load is applied. The vertical displacements of points *A *and *B *are listed in Table 2. It should be added that both points' displacements are equal to 0.3558 (Cen et al., 2009). This analysis investigates the efficiency of the elements for shear force displacements in the irregular mesh. The performances of the elements meshed irregularly in shear deflections are investigated by this test. Obviously, the error of the SST10 element's answers is approximately *2%.*

**7.1.3 Thick curved beam**

Figure 5 demonstrates the thick carved beam which is analyzed in this analysis. This structure is subjected to the shear force loading at the free end. The magnitude of the load is 600. The elasticity modulus, Poisson's ratio and the beam's thickness are 1000, 0 and 1, respectively. As before, consistence unties are utilized. By employing the suggested triangular element, this structure is meshed, as shown in Figure 5.2. It is worth emphasizing that other researchers have analyzed this beam by using quadrilateral elements. They have utilized four elements, as shown in Figure 5.1. The vertical displacement of point *A *is presented in Table 3. It should be added that the accurate displacement of the aforesaid point is 90.1 (Cen et al., 2007). The findings prove that SST10 element performs much better than well-known quadrilateral elements.

**7.1.4 MacNeal's thin cantilever beam**

This test assesses the loss of accuracy produced by utilizing elements shaped as parallelogram or trapezoid. Figure 6.1 demonstrates the thin cantilever beam with three different meshes, including rectangular, parallelogram-shaped and trapezoidal meshes. This test is proposed by MacNeal as a benchmark to evaluate the sensitivity of quadrilateral elements to mesh distortion (Macneal and Harder, 1985). Herein, this numerical analysis is deployed for SST10 triangular element for meshes shown in Figure 6.2. Three above-mentioned meshes of Macneal's beam, and the corresponding triangular meshes are demonstrated in Figure 6. Note that a big aspect ratio and distortion in parallelogram-shaped, trapezoidal and triangular meshes makes the test difficult and reliable for evaluating the efficiency of the elements. It should be added that the elasticity modulus, Poisson's ratio, and the beam's thickness are 10000000, 0.3 and 0.1, respectively. Two loading cases are applied, pure bending induced by bending moment at the free end and bending caused by unit shear force at the tip. It should be added that the exact displacements of the tip are respectively 0.0054 and 0.1081, under these load cases (Cen et al., 2009).

For mesh *A, *Table 4 indicates that SST10 element performs well when the structure is subjected to the aforesaid load cases. Furthermore, SST10 has low sensitivity to distortion of both meshes, namely *B *and *C. *On the other hand, other researchers' elements are generally sensitive to the distortion of parallelogram-shaped and trapezoidal meshes. For the trapezoidal mesh, the error of response increases extensively under the aforementioned load cases.

**7.1.5 Cantilever beam with distortion parameter**

As it is shown in Figure 7.1, a mesh which includes two quadrilateral elements is used to analyze this structure. The shape of the two elements varies with a variety of the distorted parameter e. For e=0, the elements become rectangular. However, with the increase of e, the mesh will be distorted more and more seriously. This test aims to evaluate the sensitivity of elements to the distortion of mesh. The elasticity modulus, Poisson's ratio and the beam's thickness are 1500, 0.25 and 1, respectively. The beam's free end is subjected to the bending moment M=2000. For various values of distortion parameters, the vertical displacement at the tip is presented in Table 5. Note that the near exact displacement at the tip is equal to 100 (Cen et al., 2009). This structure is analyzed by employing the proposed triangular element and the mesh shown in Figure 7.2. The results demonstrated in Table 5 prove that SST10 element is insensitive to the variations of e.

**7.1.6 Higher order patch test**

To carry out this test, a straight simple beam is considered. The length and width of the beam are 10 and 1, respectively. It should be added that the structure is subjected to pure bending. Both regular and distorted meshes are used in this numerical analysis. The meshes corresponding to the triangular and quadrilateral elements are shown in Figures 8.2 and 8.3. It should be mentioned that displacements in *x *and *y *directions are denoted by *u *and *v, *respectively. The maximum values of the displacements are presented in Table 6. Note that the exact response based on beam theory is obtained (Choi et al., 2006). It is obvious that SST10 element can reach to the good results in both regular and distorted meshes.

**7.1.7 Cook's skew beam**

This structure was employed to investigate the performance of the quadrilateral elements (Cook et al., 1989). The beam is shown in Figure 9.1. In this problem, shear displacements govern the behavior, and the distorted quadrilateral elements are deployed for meshing. It should be added that one end of this beam is clamped, and the other one is subjected to a uniformly distributed shear load *P=1. *The elasticity modulus, Poisson's ratio and the beam's thickness are 1, ⅓ and 1, respectively. Four types of mesh, including 2×2, 4×4, 8×8 and 16×16, are deployed to analyze the aforementioned beam. Figure 9.2 demonstrates the case in which 2×2 mesh and SST10 element are used for analysis of Cook's skew beam. Moreover, the vertical displacement of point *C *is indicated in Table 7. For comparison, the results of other researchers' well-known elements are presented. It is important to note that the response of GT9M8 element, for a 64x64 mesh, is considered as the near exact answer (Long and Xu, 1994).

**7.2 Nonlinear analyses**

To evaluate the performance of proposed formulation, three problems are solved in the following sections. These tests concentrate on the geometrical nonlinear behavior of structures.

**7.2.1 L-shaped plate strip**

Figure 10 shows the geometry and loading of the flat L-shaped plate strip. One end of this plate is clamped, and a uniformly distributed horizontal load is applied to the other end. The elasticity modulus and Poisson's ratio are 30000000 and 0.3, respectively. To assess the sensitivity to geometric distortion, the coarse mesh shown in Figure 10.1 is employed to analyze this structure. It should be added that the fine mesh consists of 304 square elements is shown in Figure 10.2. The results of using QM6 and Qnew elements are available for these meshing cases (Battini, 2008b). The load displacement plots for SS10, QM6 and Qnew elements are illustrated in Figure 11. It is important to mention that the horizontal displacement of the free end is considered in these plots. With the fine mesh, the three elements give exactly the same results. According to Figure 11, the SST10 element yields a good response even for the poor meshes.

**7.2.2 Cantilever beam subjected to shear force at its free end**

Figure 12 shows the geometry and loading of the cantilever beam with square cross section (Khosravi et al., 2007). The structural thickness, elasticity modulus and Poisson's ratio are 0.1, 20600000 and 0.3, respectively. Four various meshes are utilized to analyze this cantilever beam, as shown in Figure 12. The load displacement plots for SS10, Allman (3m) and LST-Ret elements are illustrated in Figures 13, 14 and 15, respectively. It is important to mention that the displacement of the free end is considered in these plots. According to the results, the responses of SST10 for the meshes, 2, 3 and 4 are very close to each others, and they have little errors. On the other hand, Figures 14 and 15 demonstrate comparisons of the SST10 results for mesh 4 with Allman (3m), and LST-Ret outcomes for meshes 1, 2, 3 and 4. These Figures show the higher efficiencies of the SST10.

**7.2.3 Thin cantilever beam under in-plane shear**

The cantilever beam with the length of 100, width of 1 and thickness of 1 is shown in Figure 16. The elasticity modulus and the Poisson ratio are equal to 1000000 and 0.3, respectively. This structure with a force at its free end is a severe test. The test was performed with regular meshes. These meshes for triangular elements are made of each rectangular mesh unit divided into two-triangular elements. The 1x100 mesh is employed for analyzing the cantilever beam. The 1x100 mesh shows that the number of the parts in the *x * and *y * directions are equal to 100 and 1, respectively. In this test, the triangular elements have not been used by other researchers. Therefore, for the sake of comparison, the responses of famous quadrilateral elements are deployed. In addition, the responses of finite-rotation for the Timoshenko beam element are utilized (Wisniewski and Turska, 2009). These are shown in Figures 17.

The SST10 solution curves for a tip of the cantilever are presented in Figures 17. The answers of some of the good elements are used in Figures 17.

**8 CONCLUSIONS**

In this paper, free formulations along with strain gradient notation were employed to introduce the optimality conditions into the finite element formulation. By satisfying the equilibrium equations, the number of nodal unknowns was decreased, and the robustness of the element was increased. Since triangular elements are suitable for modeling the structures with complex shape, the suggested scheme is quite useful. By utilizing several linear numerical tests, the high accuracy, rapid convergence, low sensitivity to geometry distortion, rotational invariance and lack of parasitic shear error were indicated for SST10 element. So far, more quadrilateral elements than triangular ones have been proposed by investigators. Consequently, the formers have been more widely used in comparison studies. Hence, the responses of the well-known quadrilateral elements were presented to demonstrate the capabilities of the proposed element in these tests. The proposed approach was without difficulty extended to the geometrically nonlinear structural analysis. In the suggested technique, the rigid motions were easily identified, and the linear formulations were simply converted into the nonlinear ones. Finally, simplicity and robustness of the proposed scheme in the geometrical nonlinear analysis of structures were clearly demonstrated numerically.

**References**

Allman, D. J., (1984). Compatible triangular element including vertex rotations for plane elasticity analysis. Computers and Structures 19(1-2):1-8. [ Links ]

Battini, J. M., (2007). A modified corotational framework for triangular shell elements. Computer Methods in Applied Mechanics and Engineering 196:1905-1914. [ Links ]

Battini, J. M., (2008a). A rotation-free corotational plane beam element for non-linear analyses. International Journal for Numerical Methods in Engineering 75:672-689. [ Links ]

Battini, J. M., (2008b). A non-linear corotational 4-node plane element. Mechanics Research Communications 35:408-413. [ Links ]

Bergan, P. G., and Felippa, C. A., (1985). A triangular membrane element with rotational degrees of freedom. Computer Methods in Applied Mechanics and Engineering 50:25-69. [ Links ]

Bergan, P. G., and Nygard, M. K., (1984). Finite Elements with Increased Freedom in Choosing Shape Functions. International Journal for Numerical Methods in Engineering 20:643-663. [ Links ]

Caylak, I., and Mahnken, R., (2011). Mixed finite element formulations with volume bubble functions for triangular elements. Computers and Structures 89:1844-1851. [ Links ]

Cen, S., Chen, X. M., and Fu, X. R., (2007). Quadrilateral membrane element family formulated by the quadrilateral area coordinate method. Computer Methods in Applied Mechanics and Engineering 196(41-44):4337-4353. [ Links ]

Cen, S., Chen, X. M., Li, C. F., and Fu, X. R., (2009). Quadrilateral membrane elements with analytical element stiffness matrices formulated by the new quadrilateral area coordinate method (QACM-II). International Journal for Numerical Methods in Engineering 77:1172-1200. [ Links ]

Chen, X. M., Cen, S., Long, Y. Q., and Yao, Z. H., (2004). Membrane elements insensitive to distortion using the quadrilateral area coordinate method. Computers and Structures 82:35-54. [ Links ]

Choi, N., Choo, Y. S., and Lee, B. C., (2006). A hybrid Trefftz plane elasticity element with drilling degrees of freedom. Computer Methods in Applied Mechanics and Engineering 195:4095-4105. [ Links ]

Choo, Y. S., Choi, N., and Lee, B. C., (2006). Quadrilateral and triangular plane elements with rotational degrees of freedom based on the hybrid Trefftz method. Finite Elements in Analysis and Design 42:1002-1008. [ Links ]

Cook, R. D., (1987). Plane hybrid element with rotational d. o. f. and adjustable stiffness. International Journal for Numerical Methods in Engineering 24(8):1499-1508. [ Links ]

Cook, R. D., Malkus, D. S., and Plesha, M. E., (1989). Concepts and Applications of Finite Element Analysis (3rd edn), Wiley (New York). [ Links ]

Dow, J. O., (1999). A Unified Approach to The Finite Element Method and Error Analysis Procedures, Academic Press. [ Links ]

Eom, J., Ko, J., and Lee, B. C., (2009). A macro plane triangle element from the individual element test. Finite Elements in Analysis and Design 45:422-430. [ Links ]

Felippa, C. A., (1994). A survey of parametrized variational principles and applications to computational mechanics. invited Chapter in Science and Perspectives in Mechanics, ed. by Nayroles, B., Etay, J., and Renouard, D., ENS Grenoble, Grenoble, France, 1-42, 1994. [ Links ] expanded version in Computer Methods in Applied Mechanics and Engineering 113:109-139. [ Links ]

Felippa, C. A., (2000). Recent advances in finite element templates. Computational Mechanics for the Twenty-First Century. Saxe-Coburn Publications, Edinburgh, pp. 71-98. [ Links ]

Felippa, C. A., (2003a). A study of optimal membrane triangles with drilling freedoms. Computer Methods in Applied Mechanics and Engineering 192:2125-2168. [ Links ]

Felippa, C. A., (2003b). A Template Tutorial, delivered in Centro Internacional de Metodos Numericos en Ingeneria, (CIMNE), Barcelona, Spain. [ Links ]

Felippa, C. A., (2006). Supernatural QUAD4: A template formulation. Computer Methods in Applied Mechanics and Engineering 195:5316-5342. [ Links ]

Felippa, C. A., and Militello, C., (1990). Variational Formulation of High Performance Finite Elements: Parametrized Variational Principles. Computers and Structures 36:1-11. [ Links ]

Felippa, C. A. and Militello, C., (1999). Construction of optimal 3-node plate bending elements by templates. Computational Mechanics 24:1-13. [ Links ]

Fredriksson, M., and Ottosen, N. S., (2004). Fast and accurate 4-node quadrilateral. International Journal for Numerical Methods in Engineering 61:1809-1834. [ Links ]

Huang, M., Zhao, Z., and Shen, C., (2010). An effective planar triangular element with drilling rotation. Finite Elements in Analysis and Design 46:1031-1036. [ Links ]

Ibrahimgovic, A., Taylor, R. L., and Wilson, E. L., (1990). A robust quadrilateral membrane finite element with rotational degrees of freedom. International Journal for Numerical Methods in Engineering 30:445-457. [ Links ]

Izzuddin, B. A., (2005). An enhanced co-rotational approach for large displacement analysis of plates. International Journal for Numerical Methods in Engineering 64:1350-1374. [ Links ]

Khosravi, P., Ganesan, R., and Sedaghati, R., (2007). Corotational non-linear analysis ofthin plates and shells using a new shell element. International Journal for Numerical Methods in Engineering 69:859-885. [ Links ]

Li, Z., and Vu-Quoc, L., (2007). An efficient co-rotational formulation for curved triangular shell element. International Journal forNumerical Methods in Engineering 72:1029-1062. [ Links ]

Long, Y. Q., and Xu, Y., (1994). Generalized conforming quadrilateral membrane element with vertex rigid rotational freedom. Computers and Structures 52(4):749-755. [ Links ]

MacNeal, R. H., (1978). Derivation of ElementStiffness Matrices by Assumed Strain Distribution. Nuclear Engineering and Design 70:3-12. [ Links ]

Macneal, R. H., and Harder, R. L., (1985). A proposed standard set of problems to test finite element accuracy. Finite Elements Analysis and Design 1(1):3-20. [ Links ]

Macneal, R. H., and Harder, R. L., (1988). A refined four-noded membrane element with rotation degrees of freedom. Computers and Structures 28(1):75-84. [ Links ]

Norachan, P., Suthasupradit, S., and Kim, K. D., (2012). A co-rotational 8-node degenerated thin-walled element with assumed natural strain and enhanced assumed strain. Finite Elements in Analysis and Design 50:70-85. [ Links ]

Pan, W., and Wheel, M., (2011). A finite-volume method for solids with a rotational degrees of freedom based on the 6-node triangle. International Journal for Numerical Methods in Engineering 27:1411-1426. [ Links ]

Pian, T. H. H., and Sumihara, K., (1984). Rational approach for assumed stress finite elements. International Journal for Numerical Methods in Engineering 20:1685-1695. [ Links ]

Piltner, R., and Taylor, R. L., (1995). A quadrilateral mixed finite element with two enhanced strain modes. International Journal for Numerical Methods in Engineering 38:1783-1808. [ Links ]

Piltner, R., and Taylor, R. L., (1997). A systematic constructions of B-bar functions for linear and nonlinear mixed enhanced finite elements for plane elasticity problems. International Journal for Numerical Methods in Engineering 44:615-639. [ Links ]

Taylor, R. L., Beresford, P. J., and Wilson, E. L., (1976). A non-conforming element for stress analysis. Int International Journal for Numerical Methods in Engineering 10:1211-1219. [ Links ]

Urthaler, Y., and Reddy, J. N., (2005). A corotational finite element formulation for the analysis of planar beams. Communications in Numerical Methods in Engineering 21:553-570. [ Links ]

Wilson, E. L., Taylor, R. L., Doherty, W. P., and Ghabussi, T., (1973). Incompatible displacement models. In: Fenven, S. T., et al. (Eds.), Numerical and Computer Methods in Structural Mechanics. Academic Press, New York, pp. 43-57. [ Links ]

Wisniewski, K., and Turska, E., (2009). Improved 4-node Hu-Washizu elements based on skew coordinates. Computers and Structures 87:407-424. [ Links ]

Yeo, S. T., and Lee, B. C., (1996). Equivalence between enhanced assumed strain method and assumed stress hybrid based on the Hellinger-Reissner principle. International Journal for Numerical Methods in Engineering 39:3083-3099. [ Links ]

Zhou, Y. X., and Sze, K. Y., (2012). A geometric nonlinear rotation-free triangle and its application to drape simulation. International Journal for Numerical Methods in Engineering 89:509-536. [ Links ]

* Corresponding author email: mrpa-jand@yahoo.com