Linear buckling analysis of thin-walled members combining the Generalised Beam Theory and the Finite Element Method

Computers and Structures 89 (2011) 1982–2000

Contents lists available at ScienceDirect

Computers and Structures journal homepage: www.elsevier.com/locate/compstruc

Linear buckling analysis of thin-walled members combining the Generalised Beam Theory and the Finite Element Method Miquel Casafont ⇑, Frederic Marimon, Magdalena Pastor, Miquel Ferrer Department of Strength of Materials and Structural Engineering, School of Industrial Engineering of Barcelona, Technical University of Catalonia, Barcelona, Spain

a r t i c l e

i n f o

Article history: Received 30 November 2009 Accepted 16 May 2011 Available online 14 June 2011 Keywords: Thin-walled members Linear buckling analysis Generalised Beam Theory Constrained Finite Strip Method

a b s t r a c t A finite element procedure to carry out linear buckling analysis of thin-walled members is developed on the basis of the existing Generalised Beam Theory (GBT) and constrained Finite Strip Method (cFSM). It allows designers to uncouple the buckling modes of a finite element model and, consequently, to calculate pure elastic buckling loads. The procedure can easily be applied to members with general boundary conditions subjected to compression or bending. The results obtained are rather accurate when compared to the values calculated via GBT and cFSM. As a consequence, it is demonstrated that linear buckling analyses can be performed with the Finite Element Method in a similar way as can be done with the existing GBT and cFSM procedures. Ó 2011 Civil-Comp Ltd and Elsevier Ltd. All rights reserved.

1. Introduction Nowadays, three different methods are usually applied to carry out linear buckling analyses of thin-walled cold-formed members: the Finite Element Method (FEM), the Finite Strip Method (FSM) and the Generalised Beam Theory (GBT). GBT presents two major advantages over the other two methods. On the one hand, it allows to take into account the distortion of the cross-section by means of beam-type finite elements [1]. Consequently, the use of GBT involves substantial reduction of the computational cost (see [2–5]). On the other hand, elastic buckling loads of pure modes (Figs. 1 and 2), or of chosen combinations of modes, can be calculated by means of GBT. This is useful under different circumstances: when pure elastic buckling loads are used to determine the strength of a member [6]; or when the contribution of each individual mode to a combined buckling mode is investigated (buckling mode identification [7,8]). Recently, the research team led by Silvestre and Camotim has developed GBTUL [9], a computer program that allows for carrying out linear buckling analyses by means of GBT. As a consequence, GBT can already be applied by practical designers of cold-formed members. Some of the concepts of GBT were translated into FSM by Ádány et al. in [10–13]. They developed a method for the calculation of pure buckling loads. The idea of the method is to force the FSM model to buckle in a given mode type by properly constraining the nodal displacements of the mesh. For this reason the method ⇑ Corresponding author. Tel.: +34 93 405 43 22; fax: +34 93 401 10 34. E-mail address: [email protected] (M. Casafont).

is called constrained Finite Strip Method (cFSM). The cFSM is implemented in CUFSM [14], an open source FSM program developed by Schafer. Pure buckling modes defined by the vector space of local, distortional or global deformation modes can be calculated with this program. Actually, CUFSM was the first program that allowed designers to carry out this sort mode decomposition analyses. The aim of the investigation presented herein is to extend the cFSM to the Finite Element Method. Since the use of FEM with shell elements is rather common in the field cold-formed steel design, the authors believe that it makes sense to try to combine this method with GBT. The investigation is developed on the basis of the work carried out by Ádány and Schafer. Pure buckling loads and modes will be obtained by constraining displacements of the finite element model. The constraining procedure consists of two steps. The first one is a GBT analysis of the cross-section of the member. Pure deformation modes, such as those shown in Figs. 1 and 2, are deduced from this analysis. Afterwards, in the second step, FEM linear buckling analyses of the member are carried out forcing the finite element mesh to deform according to the calculated individual GBT deformation modes. The procedure presented was already applied in a previous investigation to calculate pure distortional elastic buckling loads of simply supported members subjected to compression [15]. As the results obtained were satisfactory, it was decided to extend it to other modes of buckling, and to other load and boundary conditions. It is worth mentioning that a procedure for the calculation of pure buckling loads and modes by means of FEM linear bucking analyses has been recently presented by Ádány et al. in [7,8]. The

0045-7949/$ - see front matter Ó 2011 Civil-Comp Ltd and Elsevier Ltd. All rights reserved. doi:10.1016/j.compstruc.2011.05.016

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

1983

Nomenclature

a b B C E Ke Kg L n

m r R RU RW Rh

angle of the elements of the cross-section with respect to the global X axis width of the elements of the cross-section GBT basic stiffness matrix GBT basic stiffness matrix Young’s modulus finite element linear stiffness matrix finite element geometric stiffness matrix member length number of elements of the cross-section Poisson’s ratio number of half-sine waves constraint matrix U–V relationship matrix W–V relationship matrix h–V relationship matrix

RUin U–V relationship matrix for internal nodes RWin W–V relationship matrix for internal nodes t thickness of the elements of the cross-section u, v, w local displacements ui, vi, wi local node displacements Ui, Vi, Wi, hi global node displacements U; V; W; h vectors of nodal displacements m m m Um constraint factors of mode m i ; V i ; W i ; hi Vu longitudinal unknown displacement vi(x) function of the longitudinal displacement/elementary longitudinal displacement function w(x) deflection function wi(x) function of the transverse displacement in the z direction/elementary transverse displacement function. x, y, z local coordinate axes X, Y, Z global coordinate axes

main difference between this procedure and the procedure proposed in the present paper is the fact that in [7,8] the finite element mesh is not forced to buckle according to pure modes. Instead, a standard FEM linear buckling analysis is carried out and, consequently, combined buckling modes are obtained. ‘‘Pure’’ modes are identified by means of a buckling identification procedure, which provides with the percentage of participation of each mode (local, distortional or global) in the combined mode that results from the FEM analysis. It is considered that ‘‘pure’’ modes are those which show a dominant participation of a pure mode (70– 80%, see [7]). In the present investigation, the commercial finite element package ANSYS is chosen to carry out the analyses. The general formulation of the constraining method has been adapted to the analysis procedures implemented in this software. This involves some additional work that is also presented in the paper, and that may be useful for the application of the method in other commercial finite element packages. An outline of the paper follows. Section 2 presents the theoretical fundamentals of GBT needed for carrying out the cross-section analyses. Section 3 shows how to constrain FEM models to buckle in pure GBT modes. Section 4 includes an explanation on the implementation of the constraining procedure in commercial

software. An illustrative example of individual mode calculation is presented in Section 5. In Section 6, the results of additional analyses on different channel sections are shown to fully verify the proposed method. Section 7 deals with buckling loads of members with general boundary conditions. Finally, the concluding remarks of Section 8 close the paper.

Fig. 1. Longitudinal displacements of pure modes of buckling.

Fig. 2. Transverse displacements of pure modes of buckling.

2. Theoretical bases of GBT This section presents basic concepts of GBT related to the derivation of pure deformation modes of unbranched cross-sections. It

1984

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

is a summary to help readers non familiar with GBT understand the subsequent sections of the paper. Refs. [10,11,16,17] are recommended to those interested in deep knowledge of the theory. Most of the contents of the section are devoted to explain the calculation procedure of global and distortional GBT deformation modes. The calculation of local deformation modes is included at the end, in Section 2.4. 2.1. Main assumptions of GBT

where w(y) is usually taken as a sinusoidal function that provides the variation of the u(x), v(x) and w(x) mid-plane displacements along the longitudinal direction. w(y) can also be approximated by polynomials (mostly cubic). In (3.b), v(x, y) depends on the derivative of w(y) because the member shear strains are considered null (assumption (2)). Substituting w(y) by a sinusoidal function, the resultant displacement expression is:

uðx; yÞ ¼ uðxÞ  sin

In the formulation of the Generalised Beam Theory, members are considered to be composed of plates, as shown in Fig. 3. GBT assumes that these plates behave according to the Kirchhof plate theory: straight lines normal to the mid-plane remain straight, inextensional and normal to the mid-plane after deformation. Consequently, ezz = cxz = cyz = 0. Furthermore, two additional simplifying assumptions are considered: 1. The membrane transverse extensions are zero:

eMxx ¼ 0:

ð1Þ

v ðx; yÞ ¼

M xy

c ¼ 0:

ð2Þ

These are the two main assumptions of GBT (see [16,17]). Fig. 3 shows the notation that will be used for the displacements of the plate mid-plane points (the notation given in [10,11] is followed): u, v, w and h are the displacements expressed in the local plate systems, the x–y–z coordinate systems; and U–V– W and h are the displacements corresponding to the global coordinate system, the X-Y-Z system. Later on, it will be seen that it is very important to make the distinction between longitudinal displacements: v and V; and transverse displacements: u, w, U, W and h. In GBT, as in FSM, displacements u(x, y), v(x, y) and w(x, y) are expressed as a product of two single-variable functions:

uðx; yÞ ¼ uðxÞ  wðyÞ;

ð3:bÞ

wðx; yÞ ¼ wðxÞ  wðyÞ;

ð3:cÞ

ð4:bÞ

r  p  y ; L

ð4:cÞ

where L is the length of the member and r the number of half-waves.One of the most important features of GBT is that transverse displacements u(x) and w(x) can be expressed in terms of the longitudinal displacements v(x). Firstly, v(x) is expressed in terms of the longitudinal node displacements vi, nþ1 X

v i ðxÞ  v i ;

ð5Þ

i¼1

where n is the number of elements of the cross-section (Fig. 4a), vi is the longitudinal displacement of node i, and vi(x) is a linear function of x. It should be observed that only nodes located in the wall longitudinal edges are considered. These nodes are called natural nodes. Additional intermediate nodes will be needed when dealing with local modes (see Section 2.4). Subsequently, what is actually done is to express the transverse displacements in terms of the longitudinal node displacements:

uðx; yÞ ¼

nþ1 X

ui ðxÞ  v i  sin

i¼1

v ðx; yÞ ¼

ð3:aÞ

v ðx; yÞ ¼ v ðxÞ  w;y ðyÞ;

ð4:aÞ

r  p  y rp ;  v ðxÞ  cos L L

wðx; yÞ ¼ wðxÞ  sin

v ðxÞ ¼

2. The membrane shear strains are zero:

r  p  y ; L

¼

nþ1 r  p  y X ¼ ui ðxÞ  /i ðyÞ; L i¼1

ð6:aÞ

nþ1  r  p  y X rp v i ðxÞ  v i  cos L L i¼1 nþ1 X

v i ðxÞ  /i;y ðyÞ;

ð6:bÞ

i¼1

wðx; yÞ ¼

nþ1 X

wi ðxÞ  v i  sin

 r  p  y

¼

L

i¼1

nþ1 X

wi ðxÞ  /i ðyÞ;

i+1

Y, V

ð6:cÞ

i¼1

i+2 n-1 n

i

b

i

n

θ

x, u i

n

X, U

n n+1

i+1

n+1

i

n+1

i z, wi

i+1

X, U

i

X, U

x, u y, v

Z, W

i θ

1

z, w

1

2

Z, W

2 3

Z, W

2

3

(a) Fig. 3. Notation for the cross-section elements, nodes and degrees of freedom.

1

1 3

2

(b)

Fig. 4. Nodes considered in the formulation of GBT: (a) natural nodes; (b) natural and intermediate nodes.

1985

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

where ui(x) and wi(x) are functions of the cross-section geometry and material properties of the member. These functions can be determined through a rather long procedure explained in [16,17]. In the present paper, only the two functions needed for the calculation of the GBT deformation modes are derived: – relationship between transverse and longitudinal nodal displacements, presented in Section 2.2; – relationship between the function of transverse displacements wi(x) and the longitudinal displacements for a particular vi(x), presented in Section 2.3. 2.2. Relationship between transverse and longitudinal node displacements The transverse-longitudinal displacement relationships are determined by following a two step procedure. In the first step, the global transverse displacements Ui and Wi of the internal nodes (nodes 2 to n in Fig. 4a), are expressed in terms of the global longitudinal displacements Vi. This relationship, which is called here UWin–V, is derived from the GBT main assumptions (1) and (2), and from special geometric compatibility conditions. In the second step, the equilibrium condition of the cross-section is applied to determine: – UWex–V relationship between Vi and the U1, Un+1, W1 and Wn+1 transverse displacements of the external nodes (nodes 1 and n + 1 in Fig. 4a); – h–V relationship between Vi and hi rotations of all nodes. 2.2.1. Derivation of the uwin–v relationship When the Generalised Beam Theory is derived, the variation of u(x) and v(x) in Eq. (5) is assumed to be linear within the elements of the cross-section [17]. For example, in [10] these displacements are expressed in terms of the nodal degrees of freedom ui and vi and a linear function of x:

   r  p  y x x  ui þ  uiþ1  sin ; uðx; yÞ ¼ 1 bi bi L

v ðx; yÞ ¼

rp  L

ð7:aÞ

   r  p  y x x  v i þ  v iþ1  cos : 1 bi bi L

ð7:bÞ

If the simplifying assumption (1) is taken into account, the following condition is obtained:

eMxx ¼

@u uiþ1  ui ¼0!  sin @x bi

 r  p  y L

¼0

v i  v iþ1

ui ¼

bi

ð12Þ

:

As vi = Vi (see Fig. 3), displacement ui can be directly related to the longitudinal global displacement:

V i  V iþ1 : bi

ui ¼

ð13Þ

Once the relationship corresponding to the ui component is known, geometric compatibility conditions provide an expression of wi in terms of vi. This can be seen in Fig. 5, where an arbitrary longitudinal vi displacement is applied at node i, considering at the same time that, for the moment, the cross-section elements are connected by means of perfect hinges. Actually, the displacement wi is first expressed in terms of the transverse displacement ui, taking into account that there must be continuity of the crosssection at nodes. Afterwards, (12) (or (13)) is used to express wi in terms of vi (or Vi). The relationship between wi and ui can be deduced from the drawing in Fig. 6, which shows an arbitrary displacement of an internal node of the section. The displacement has been drawn in such a way that the wi–ui relationship is easily derived. However, the final equation is valid for any arbitrary displacement. This relationship is obtained with the help of the two triangles highlighted in Fig. 6b:

wi ¼

ui1 ui  ; sinðDai Þ tanðDai Þ

ð14Þ

where:

Dai ¼ ai  ai1 :

ð15Þ

At this point, the relationships between the transverse ui and wi displacements and the longitudinal vi (or Vi) displacements have already been obtained (Eqs. (12) and (14)). 2.2.2. Derivation of the UWin–V relationship It will be shown in Section 3 that the GBT constraints will be applied to the finite element mesh in the global coordinate system. As a consequence, equation (14) has to be translated into global coordinates. This can be done by means of the following equations (see Fig. 6a):

U i ¼ ui  cosðai Þ þ wi  sinðai Þ   ui1 ui ¼ ui  cosðai Þ þ   sinðai Þ; sinðDai Þ tanðDai Þ

ð16:aÞ

ð8Þ

Since (22) has to be fulfilled in all cross-sections of the member, then:

uiþ1 ¼ ui

i+1

ð9Þ bi

Therefore, the u(x, y) displacement within an element of the cross-section can be expressed as:

r  p  y uðx; yÞ ¼ uðyÞ ¼ ui  sin : L

u i =vi /bi

vi

u i-1 =vi /bi-1

ui i

r  p  y @u @ v rp þ ¼0!  ui  cos @y @x L L  r  p  y r  p v iþ1  v i ¼ 0;  cos  þ L L bi

u i-1 wi-1

b i-1

1

ð11Þ

2 dy

and, consequently, the relation between longitudinal transverse ui displacements is obtained:

n+1

wi

ð10Þ

Furthermore, if the simplifying assumption (2) is considered, then:

cMxy ¼

n

i

vi

and Fig. 5. Geometric determination of the wi displacement.

1986

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

X

X

i+1

i+1 Z

Z αi

ui-1

wi-1

wi

ui

u i-1 Ui

Δα i

wi-1

wi Wi

cosðai Þ V i1  V i cosðai1 Þ V i  V iþ1 þ   sinðDai Þ sinðDai Þ bi1 bi   cosðai Þ cosðai Þ cosðai1 Þ  Vi ¼  V i1 þ þ sinðDai Þ  bi1 sinðDai Þ  bi1 sinðDai Þ  bi cosðai1 Þ   V iþ1 : ð20:bÞ sinðDai Þ  bi

Wi ¼ 

i

i

Δαi

The above equations can be expressed in matrix form:

Δα i

ui

ui-1

Ui ¼

Δα i

(

sinðai Þ sinðDai Þ  bi1



sinðai Þ sinðai1 Þ   sinðDai Þ  bi1 sinðDai Þ  bi



8 9 ) > < V i1 > = ;  Vi > > : ; Vi þ1 ð21:aÞ

α i-1 i-1

sinðai1 Þ sinðDai Þ  bi

8 9 (   ) > < V i1 > = cosðai Þ cosðai Þ cosðai1 Þ cosðai1 Þ  Wi ¼  þ Vi ;  > > sinðDai Þ  bi1 sinðDai Þ  bi1 sinðDai Þ  bi sinðDai Þ  bi : ; V iþ1

i-1

Fig. 6. Geometric relationship between transverse displacements.

ð21:bÞ

W i ¼ ui  sinðai Þ  wi  cosðai Þ   ui1 ui ¼ ui  sinðai Þ    cosðai Þ: sinðDai Þ tanðDai Þ

ð16:bÞ

and they can be extended to all internal nodes of the cross-section:

     3 2 sinða2 Þ a2 Þ a1 Þ sinða1 Þ 8 8 9 9  sinðsinð  sinðsinð 0 0 0 0 0 Da2 Þb1 Da2 Þb2 sinðDa2 Þb1 sinðDa2 Þb2 V1 > U 7 > 6 > > 2 > > > >       7 > > > > > > > 6 > sinða3 Þ a3 Þ a2 Þ sinða2 Þ > > 6 U3 > V2 > > > > >  sinðsinð 0 0 0 07  sinðsinð 0 > > > > 7 6 D a Þb D a Þb D a Þb D a Þb sinð sinð 3 2 3 2 3 3 3 3 < < = 6 = 7 .. . 7 6         . ¼ . > 6 7 > . > > > > > > > > 6        7 > > Vn > > 7 > > U n1 > > 6 > > > > > 7 > 6 > > > > : ; 4 ; 5 : 0 0 0 0 x x x 0 Un V nþ1 0 0 0 0 0 x x x

ð22:aÞ

     3 2 cosða2 Þ cosða2 Þ cosða1 Þ cosða1 Þ 8 8 9 9  sinð  sinð 0 0 0 0 0 sinðDa2 Þb1 sinðDa2 Þb2 Da2 Þb1 Da2 Þb2 V1 > W 7 > 6 > > 2 > > > >     7 6 > > > > > > > > cosða3 Þ cosða3 Þ cosða2 Þ cosða2 Þ 7 > 6 > W3 > V2 > > > > > Þ   0 ðsinð 0 0 0 0 > > > > 7 6 D a Þb D a Þb D a Þb D a Þb sinð sinð sinð 3 2 3 2 3 3 3 3 < < = = 7 6 .. .. 7 6         ¼  . . 7 6 > > > > > > > > 6 > >        7 > > Vn > > > 7 > 6 > > > W n1 > > > > 7 > 6 > > > > : ; ; 5 : 4 0 0 0 0 x x x 0 Wn V nþ1 0 0 0 0 0 x x x

ð22:bÞ

After some trigonometric operations performed taking into account that:

or,

sinðDai Þ ¼ sinðai Þ  cosðai1 Þ  cosðai Þ  sinðai1 Þ;

ð17Þ

U in ¼ RUin  V;

ð23:aÞ

cosðDai Þ ¼ cosðai Þ  cosðai1 Þ þ sinðai Þ  sinðai1 Þ;

ð18Þ W in ¼ RWin  V:

ð23:bÞ

the following relations are obtained:

Ui ¼

sinðai Þ sinðai1 Þ  ui1   ui ; sinðDai Þ sinðDai Þ

Wi ¼ 

cosðai Þ cosðai1 Þ  ui1 þ  ui : sinðDai Þ sinðDai Þ

ð19:aÞ

ð19:bÞ

Introducing the relationship between the transverse ui displacements and the longitudinal Vi displacements (equation (13)) in the above equations, the final relation is achieved:

sinðai Þ V i1  V i sinðai1 Þ V i  V iþ1    sinðDai Þ sinðDai Þ bi1 bi   sinðai Þ sinðai Þ sinðai1 Þ  Vi ¼  V i1  þ sinðDai Þ  bi1 sinðDai Þ  bi1 sinðDai Þ  bi sinðai1 Þ þ  V iþ1 ð20:aÞ sinðDai Þ  bi

Ui ¼

2.2.3. The equilibrium of the cross-section In the present section, the UWex–V and h–V relationships are derived. They are determined through the relation between the internal displacements Uin and Win and the other transverse displacements of the cross-section (from now on in-exh relationship). Once in-exh is known, equations (23) can be used to obtain the final relation with the global longitudinal degrees of freedom. The in-exh relationship is determined by means of the equivalent 2D beam model shown in Fig. 7. The analogy between the cross-section and a beam system is applied, assuming that the cross-section deforms so that this 2D system is in equilibrium. Now, the connection between the cross-section elements of the model is perfectly rigid, and the internal nodes are assumed to be supported by means of external hinges where displacements Uin and Win have been imposed.

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

The following steps should be applied to obtain the in-exh relation (similar to the standard 2D beam finite element procedure):

5. The system of equations of the global equilibrium is reordered to separate the ‘‘known’’ displacements from the ‘‘unknown’’ displacements:

1. For each beam element (Fig. 7), the stiffness matrix is calculated in the local coordinate system:

2

EAi bi

6 6 0 6 6 6 6 0 6 K i ¼ 6 EA i 6 6 bi 6 6 0 6 4 0

0

0

EAi bi

0

12EIi b3i

6EIi b2i

0

i  12EI b3

6EIi b2i

4EIi bi

0

0

EAi bi

0

0

0 i

i  6EI b2 i

i  12EI b3

i  6EI b2

0

12EIi b3i

6EIi b2i

2EIi bi

0

i  6EI b2

i

i

i

3

7 7 7 7 2EIi 7 7 bi 7 7; 0 7 7 7 i 7  6EI b2i 7 5 6EIi b2i

" ð24Þ

K kk K uk

EAi bi

where:

ti ð1  t2 Þ

ð25Þ

t 3i : 12  ð1  m2 Þ

ð26Þ

Ai ¼

Ii ¼

2. The element stiffness matrices are expressed in the X–Z global coordinate system: T

Ki ¼ T  Ki  T

ð27Þ

where T is the usual coordinate transformation matrix. 3. Assembly of the global stiffness matrix K from K i . 4. Global equilibrium:

K d¼K 

8 > > > > > > > > > > > > > > > > > > > <

U1 W1 h1 U2 W2

8 9 0 > > > > > > > > > 0 > > > > > > > > > > > 0 > > > > > > > > > F > > U2 > > > > > > < F W2 > =

9 > > > > > > > > > > > > > > > > > > > =

¼ h2 > > > > > > > > > > > .. > > > > > > > > > > . > > > > > > > > > > > > > > > U > > > nþ1 > > > > > > > > > > > > W > > > nþ1 > > : : ; > hnþ1

0 .. .

> > > > > > > > > > > 0 > > > > > 0 > > > ; 0

8 9 8 9 U2 > F U2 > > > > > > > > > > > > > .. > .. > > > > > > > > > > > > > . . > > > > > > > > > > > > > > > > > > > Un > F Un > > > > > > > > > > > > > > > > > > > > > > W2 > F W2 > > > > > > > > > > > > > > > > > . . > > > > . . > > > > > > > # > . . > > > > < < = = K ku W F n Wn  ¼ ; > > > > > K uu > > > > >  >  > > > > > > > > > > > > > > > > > > > U1 > 0 > > > > > > > > > > > > > > > > > > > > > 0 U > > > > nþ1 > > > > > > > > > > > > > > > > 0 W > > > > 1 > > > > > > > > > > > > > > > > 0 W > > > > nþ1 > > > > > > > > : ; : ; h 0

ð29Þ

or in more compact form according to the notation given by Ádány in [18]:

"

K kk K uk

8 9 8 9 # > d k > > qk > < = < = K ku   ¼  ; : > ; : > ; > K uu > 0 du

ð30Þ

where dk is the vector of ‘‘known’’ displacements, du is the vector of ‘‘unknown’’ displacements, qk may be seen as the reaction forces provoked by the imposed ‘‘known’’ displacements, 0 is a null vector, K kk ; K uk ; K ku and K uu are constructed from the global stiffness matrix K of Eq. (28). 6. The in-exh relationship is calculated from (30):

du ¼ K 1 uu  K uk  dk ;

ð31Þ

8 9 U2 > > 8 9 > > > > > .. > U1 > > > > > > > > > > > > . > > > > > > > > > > < U nþ1 > = = n 1 W1 ¼ K uu  K uk  ; > > > W2 > > > > > > > > > W > > > > nþ1 > > > .. > > > > > : > ; > > > > h > . > > : ; Wn

ð32Þ

or:

¼ q:

ti

i+1 n n+1

bi

1987

ð28Þ

and in terms of the longitudinal displacements:

8 9 8 9 U1 > > > > V1 > > > 2 3 > > > > > > > > > U > > > nþ1 R Uin < V2 > = < = 6 7 1 W1 : ¼ K uu  K uk  4    5  . .. > > > > > > > > > > W nþ1 > > > > > > RWin : ; > > > > : ; V nþ1 h

ð33Þ

The above equation is finally expressed as follows:

i

8 9 2 3 8 9 RU1 U1 > > > > V1 > > > > > > > 7 > 6 > > > > R U > > Unþ1 7 > < nþ1 > = 6 7 < V2 = 6 7 6 W1 ; ¼ 6 RW1 7  . > > > .. > > > > 6 7 > > > > > W nþ1 > > > 4 RWnþ1 5 > : > ; > > > : ; V nþ1 h Rh 1 2 3 Fig. 7. Cross-section beam model.

ð34Þ

where RU1 , RUnþ1 , RW1 , RWnþ1 and Rh , are vectors and matrices derived from equation (33). Eq. (34) is the last step of the procedure for the determination of the transverse–longitudinal relationships. The present section will be closed with a compact equation that provides the final relation-

1988

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

ship between all transverse displacements and longitudinal displacements. This equation can be obtained by re-ordering the vectors and matrices shown in Eqs. (23) and (34)

3 2 8 9 RU1 U1 > > 7 6 > > > > > > 6 RUin 7 > 7 8 6 > U in > > 9 2 9 > 3 8 8 9 > 7 6 > > > V1 > V1 > 6 RUnþ1 7 > > > > U nþ1 > RU U > > > > > > > > > > > > > > > 7 > 6 > > > > > > > > 7 6 < V = < V = 6 I 7 < V2 = 6 7 < V2 = I 7 6 7 6  ¼ ¼ ¼6 7 .. > 6 .. >; 7 > > W1 > 6 RW1 7 > W> > > > > > > > 4 RW 5 > . . > > > > > > > > > 7 6 : ; > > > > > 6 : : ; ; > 7 > > W in > > h > > 7 6 R V V R Win nþ1 nþ1 > > h > > 7 6 > > > W nþ1 > 7 6 > > > > : ; 4 RWnþ1 5 h Rh ð35Þ where I is the (n + 1) identity matrix. 2.3. GBT cross-section analysis The derivation of the relationship between longitudinal and transverse displacements shown in the previous section can actually be considered the first part of the cross-section analysis. The method to obtain this relationship is similar in both GBT and cFSM. On the contrary, different procedures are applied to determine the individual deformation modes, which is the second part of the cross-section analysis. Since GBT individual modes are used in the present investigation, the GBT procedure for this second part of the analysis is summarised here (a detailed explanation on the cFSM procedure can be found in [11]). In GBT, the individual modes are determined through the solution of a particular eigenvalue problem where two matrices are involved, B and C. This section is devoted to show how these matrices are defined. Their calculation is carried out from two displacement functions that are presented in the next two subsections. 2.3.1. Elementary longitudinal displacement functions The elementary longitudinal displacement function, vi(x), is a linear function of the transverse coordinate x that has a unit value at node i and zero value at all other nodes (see an example in Fig. 8). It can be expressed mathematically with the following equations

8 x 6 xi1 ; > > > > < xi1 < x 6 xi ;

v i ðxÞ ¼ >

xi 6 x < xiþ1 ; > > > : x P xiþ1 ;

where xi ¼

v i ðxÞ ¼ 0; v i ðxÞ ¼ b 1  x  bx v i ðxÞ ¼  b1  x þ xb v i ðxÞ ¼ 0;

i1

i1

i1 iþ1

i

i

w0i ðxÞ ¼ w0i  u1 ;

ð36Þ

;

0 i+1

0 0

bj

d0i ðxÞ ¼

bi  h0i  u2 ; 2

xi 6 x 6 xiþ1

ð38Þ

where h0i is also determined from the local transverse displacements of nodes i and i + 1; and function u2 is

u2 ¼ 2  n  1

ð39Þ

x  xi bi

ð40Þ

n¼

Before the third part of the local transverse displacements is determined, the node bending moments, mi and mi+1, provoked by the unit longitudinal displacement of node i should be calculated. The beam finite element model of the cross-section presented in Section 2.2.3 can be used to obtain these bending moments. Afterwards, the equation to be applied for the deflections can be derived from classical beam theory 2

wdi ðxÞ ¼ 

bi  ðmi  u3 þ miþ1  u4 Þ; 3  Ki

xi 6 s 6 xiþ1 ;

ð41Þ

3 2

1 2

ð42Þ

u4 ¼   n þ  n3 ;

ð43Þ

u3 ¼ n þ  n2   n3 ;

Ki ¼

1 2

E  t3 : 12  ð1  m2 Þ

ð44Þ

The total elementary local flexural function is obtained by summation of (37), (38), and (41)

i 0

ð37Þ

where u1 = 1 and w0i is determined from the local transverse displacements of nodes i and i + 1 caused by the unit longitudinal displacement applied at node i. Fig. 5 show how these transverse node displacements can be calculated. Eqs. (12) and (14) can be used for the transverse displacement of node i, while simple additional geometrical calculations should be applied to determine the transverse displacements of node i + 1 (and i  1). The equation for the part of the transverse displacement due to the rigid body rotation is

1 2

j¼1 bj .

0

xi 6 x 6 xiþ1 ;

where

;

Pi1 1

2.3.2. Elementary local flexural functions The elementary local flexural functions describe the deformation experienced by the elements of the cross-section when the unit longitudinal displacement presented in the previous section is applied to one node (Fig. 8). Fig. 9, shows this deformation, which can be decomposed in three parts: rigid body transverse displacement (Fig. 9 b), rigid body rotation (Fig. 9c) and bending deflections (Fig. 9d). The rigid body transverse displacement is expressed in terms of the displacement of the element midpoint

i 0

2

0 0

wi ðxÞ ¼ w0i  u1 þ

0

bi b  h0i  u2  i  ðmi  u3 þ miþ1  u4 Þ: 2 3  Ki

ð45Þ

2.3.3. Calculation of the GBT deformation modes The GBT deformation modes result from the solution of the above mentioned eigenvalue problem

(a)

(b)

Fig. 8. Unit displacement at node i: (a) elementary longitudinal displacement function of node i. (b) Transverse deflections associated with the unit displacement of node i.

ðB  k  ECÞ  v ¼ 0:

ð46Þ

The components of matrices B and C (Bik and Cik) are determined from functions (36) and (45). First all elementary functions should

1989

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

i+1 ui

i+1

wi+1

ui

i

i

θ 0i

i

i+1

w0i (x)

w0i (x)

w0i (x)

i

w(x) i

ui

δ0i(x)

δ0i(x)

wdi (x)

wi

i

i

(a)

i

(b)

i

(c)

(d)

Fig. 9. Transverse deflection of element i caused by a longitudinal unit displacement of node i: (a) elementary local flexural functions, (b) rigid body displacement, (c) rigid body rotation, and d) transverse deflection.

be obtained considering a unit displacement at each node of the cross-section. Afterwards, the following equations should be applied

Bik ¼

n X j¼1

C ik ¼

n X

E  12  ð1  m2 Þ

E

Z

Z

xjþ1

! t 3  wi;xx  wk;xx  dx

xjþ1

t  ui  uk  dx þ

xj

j¼1

ð47Þ

xj

E  12  ð1  m2 Þ

Z

xjþ1

! t3  wi  wk  dx

xj

ð48Þ d2 wj

d2 wi dx2

where wi;xx ¼ and wj;xx ¼ dx2 . The important point is that the eigenvectors v ð¼ VÞ resulting from the simultaneous diagonalization of C and B have physical meaning: they are the longitudinal displacements of the GBT deformation modes. n + 1 eigenvectors are obtained, n  3 of which contain the longitudinal displacements of the distortional modes (Fig. 1). The other four eigenvectors correspond to global buckling modes. The calculation of these global eigenvectors involves more work, as it is explained in [16,17]. They can also be easily calculated from the coordinates of the natural nodes [11]. It is noted that only the longitudinal displacement components are obtained. Eq. (35) is applied to determine the transverse components. For example, the complete vector of nodal displacements of a specific GBT deformation mode m will be calculated as follows:

8 > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > <

Um 1 Um 2 .. . Um nþ1 Vm 1 Vm 2 ...

9 > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > =

3 8 m 9 2 8 m9 > V > > RU U > > 1m > > > > > > > > 7 > 6 < V2 > < Vm > = = 7 6 Vm I nþ1 m 7 ¼ ¼RV ¼6 m . 7 6 m> > > > > . W1 > W > > > > 4 RW 5 > . > > > > > > : m > > ; > > m > > > > > : ; > > m W h > > 2 V Rh > > nþ1 > > > > > > .. > > > . > > > > > > > > m > > > > > W > > nþ1 > > > m > > > > > h > > 1 > > > > m > > > > h > > 2 > > > > > > > .. > > > . > > > > > > > : m > ; hnþ1

2.4. GBT cross-sectional analysis with intermediate nodes The inclusion of intermediate nodes in the GBT cross-sectional analysis (see Fig. 4b) leads to two major improvements in the calculation procedure of deformation modes: (1) local modes can be determined, and (2) displacements at intermediate points of the cross-section elements are obtained for distortional and global modes. As it will be seen in Section 6, these intermediate displacements are needed to properly constraint a finite element mesh. However, the consideration of the intermediate nodes makes the calculation procedure presented in the previous sections more complex. So far the calculations of the GBT cross-section analysis have been carried out in the vector space defined by the unit longitudinal displacement functions presented in Section 2.3 (Fig. 8). When intermediate nodes are taken into account, the main conceptual change is that this vector space is expanded to include unit transverse displacement functions such as those shown in Fig. 10. It should be observed that these functions do not involve longitudinal displacements. The following points summarise how the general procedure changes to take into account intermediate nodes: 1. As pointed out above, elementary transverse displacement functions wi(x) are defined considering unit local displacements wi. 2. The 2D beam finite element model presented in Section 2.2.3 should be applied again to determine the nodal bending moments of equation (45), that are used in the calculation of wi(x). In this analysis, all transverse nodal displacements, but wi, are fixed to zero for both natural and intermediate nodes. 3. The elementary longitudinal and transverse displacement functions defined for distortional and global modes are valid. However, they should be divided according to the location of intermediate nodes. 4. The number of resultant eigenvectors, or GBT deformation modes, is equal to the total number of nodes plus two. Two additional modes are obtained because the nodes of the external edges (1 and n + 1 in Fig. 4) are considered natural and intermediate at the same time.

ð49Þ

m m m In the present investigation, components U m i ; V i ; W i and hi are called constraint factors of mode m. They are defined independent from the length of the member.

3. Constraining a finite element mesh Elastic buckling loads and modes are determined in GBT by solving an eigenvalue problem that is derived from the second order equation of equilibrium [19]. The equation of equilibrium, and consequently the eigenvalue problem, can be expressed in the vectorial space defined by the GBT deformation modes presented in the previous section. Instead of using displacement degrees of freedom, GBT deals with deformation modes. This makes very easy to calculate buckling loads of pure modes or chosen combinations of modes. The designer only has to choose which deformation mode or combination of modes are to be considered in the analysis,

1990

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

0

0

0

i+1

1

0

i

of a function depending on y (see Eqs. (4)). For the moment, members are considered simply supported and, consequently, a simple sinusoidal function is sufficient to constraint the nodal displacements

i+1 0 i

0

0 i-1

wi (x)

i-1

UðyÞ ¼ U  wðyÞ ¼ U  sin VðyÞ ¼ V  w;y ðyÞ ¼

0

r  p  y ; L

r  p  y rp ;  V  cos L L

WðyÞ ¼ W  wðyÞ ¼ W  sin 0 0

0 0

hðyÞ ¼ h  wðyÞ ¼ h  sin 0

0 0

0

(b)

Fig. 10. Transverse deflections caused by a transverse unit displacement of intermediate node i: (a) unit local displacement wi and (b) elementary transverse displacement function.

and the eigenvalue problem is constructed and solved for the selected modes. A completely different approach is used to calculate buckling loads in the constrained FSM. As mentioned earlier, the idea is to force the FSM models to buckle in a given mode type by properly constraining the nodal displacements of the mesh. This section is devoted to explain the FSM constraining procedure and its extension to FEM. When linear FSM buckling analyses are performed, an eigenvalue problem is also solved:

ðK e  k  K g Þ  D ¼ 0;

ð50Þ

where K e is the elastic stiffness matrix, K g is the geometric stiffness matrix, the eigenvalues k are load factors, from which the elastic buckling loads are determined, and the eigenvectors D provide the modes of buckling.If pure elastic buckling loads are to be calculated, the eigenvalue problem transforms into (see [18])

  RT  K e  R  k  RT  K g  R  DM ¼ 0;

ð51Þ

where R is the constraint matrix. For example, if the calculation is constrained to global modes, the constraint matrix is

V g2

V g3

V g4 :

ð52Þ

V gi are the longitudinal displacements of the orthogonal base vectors of the global sub-space. In the present investigation, global GBT deformation modes are used. If the calculation is constrained to an individual buckling mode, R should contain only one deformation mode. For instance, if the problem is constrained to the symmetric distortional mode, the constraint matrix will be

Rsd ¼ R  V sd

r  p  y ; L

r  p  y : L

ð54:bÞ ð54:cÞ ð54:dÞ

If Eqs. (54) are introduced in Eq. (49), the original relationship becomes valid for FEM models

(a)

RG ¼ R  HG ¼ R  ½ V g1

ð54:aÞ

ð53Þ

In this paper, the calculation of pure buckling loads and modes is carried out with the Finite Element Method applying an approach similar to the constrained Finite Strip Method. The main conceptual change introduced is that cross-section constraints, such as (52) and (53), should be extended in the longitudinal direction. This can be done by taking into account that in GBT the nodal displacements are actually expressed in terms

3 8 9 RU  sinðrpL yÞ > V1 > > 7 > 6 > > > > > 7 > 6 > > > 6 rp  I  cosðrpyÞ 7 > V < = 2 7 6 L L 7 6  ¼ RðyÞ  V; ¼6 7 .. > > > > > 6 R  sinðrpyÞ 7 > WðyÞ > > > > > . > W > > > > 7 6 L > > > > > > > 5 > > : > ; 4 : ; hðyÞ rpy V nþ1 Rh  sinð L Þ 8 9 UðyÞ > > > > > > > > > > > < VðyÞ > =

2

ð55Þ

Each cross-section along the length of the member will have its own relationship. As a consequence, the vector of constraint factors will also depend on y. For example, for mode m

8 9 rpy Um > > 1  sinð L Þ > > > > > > m rpy > > > >  sinð Þ U > > 2 L > > > > > > > > . > > . > > > > . > > > > > > m > rpy > > U nþ1  sinð L Þ > > > > > > > > m rpy > rp > > > >  V  cosð Þ > > 1 L L > > > > m > > r p y rp > >  V  cosð Þ > > 2 L L > > > > > > > > .. > > > > > > . > > > > > > > < rp  V m  cosðrpyÞ > =

8 9 > U m ðyÞ > > > > > > < V m ðyÞ > = nþ1 L L ¼ m rpy m > W ðyÞ > > W  sinð Þ > > > 1 L > > > > > : m ; > > >  sinðrpL yÞ Wm > h ðyÞ 2 > > > > > .. > > > . > > > m > > W nþ1  sinðrpL yÞ > > > > rpy > > hm 1  sinð L Þ > > > m > > h2  sinðrpL yÞ > > > > > .. > > > . > > > : m hnþ1  sinðrpL yÞ

> > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > ;

ð56Þ

A new constraint matrix R is defined from (56). It can be observed that this vector only contains GBT degrees of freedom, i.e., those degrees of freedom which are involved in the derivation of the GBT deformation modes (and in the constraining procedure of cFSM): U (or u), V, W (or w) and h. There are two degrees of freedom per node, two rotations (hx and hz), that are not constrained. Finally, it is emphasised that finite element meshes will be constrained to GBT deformation modes. This is another conceptual difference between the procedure presented here and the constrained Finite Strip Method. In cFSM, the mesh is constrained to pure sub-spaces (local, distortional or global) defined by means of base vectors that can be different from the GBT modes [11–13].

1991

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

4. Practical application of the constraining procedure For practical reasons, in the present investigation it was decided to make the constraints partially independent from the finite elem m m ment mesh. The constraint factors (U m i ; V i ; W i and hi ), whose value depends on the node location in the cross-section, are defined by means of a program in MATLAB before the finite element mesh is generated. Afterwards, the mesh generator of ANSYS is forced to create nodes in points for which the constraint factors have been previously defined, as it is shown in Fig. 11. This procedure is easier than to calculate the constraint factors from the location of the nodes of a previously defined mesh. A small change should be introduced in the transformation matrix R to take into account that, as shown in Fig. 11, not all degrees of freedom of the mesh are constrained. The new R is composed of four sub-matrices: (1) in the upper left part, there is an identity sub-matrix that will multiply the components of the stiffness matrices of (51) corresponding to the non-constrained degrees of freedom; (2) in the lower right part, the original R matrix for the constrained degrees of freedom; and (3) two zero submatrices in the upper right and lower left parts. The stiffness matrices of the eigenvalue (51) have also to be modified, reordering their components according to the order of the constraint matrix, as explained in [15]. However, in this investigation matrix R is not used because ANSYS does not allow for modifying the stiffness matrices according to Eq. (51). In ANSYS, the mesh is forced to buckle in a specific deformation mode by introducing constraint equations between degrees of freedom. Once the constraint equations are set, the program solves the constrained problem applying the Lagrange multipliers method [20]. This method is always adding equations to the global system, one per constraint equation, unlike the simple condensation procedure, by which the condensed degrees of freedom are eliminated. Therefore, the larger the number of constraint equations, the larger the calculation time. The next subsection shows how these constraint equations are derived for the case in which only one GBT deformation mode is considered. Later on, in Section 4.2, it will be shown the derivation of the constraint equations for combined GBT modes.

strained degree of freedom. A longitudinal degree of freedom, Vu in Fig. 11.a, is chosen to be the unknown degree of freedom. Eq. (56) provides the relationship between the unknown degree of freedom and the other constrained displacements. If the member is forced buckle according to mode m, the longitudinal components of the constrained nodes should accomplish the relationship given by (56), i.e., the displacement vector of these nodes has to be proportional to (56) by a factor of C

VðyÞ ¼

Vu

9 > > > > > > > > > > > > > =

rpy > rp >  Vm > u  cosð L Þ > L > > > > > > > > > > .. > > > > > > . > > > > > > : rp m rpy ;  V nþ1  cosð L Þ L

 C:

ð57Þ

Factor C can be expressed in terms of the unknown constrained degree of freedom by considering that

r  p  y  rp m u  V u  cos  C; L L

Vu ¼

ð58Þ

where yu is the longitudinal coordinate of the node to which the unknown degree of freedom belongs. An expression for C is obtained from (58)

C ¼ rp L

Vu rpyu  Vm  u cosð L Þ

ð59Þ

Finally, the relationship between Vu and the other constrained nodes is obtained substituting (59) in (57)

VðyÞ ¼

4.1. Constraint equations for individual GBT deformation modes Constraint equations are derived only for the longitudinal degrees of freedom to easily illustrate the procedure. The first step is to choose the degree of freedom to which all constrained degrees will be linked. This degree of freedom is called the unknown con-

8 rp m  V 1  cosðrpL yÞ > L > > > rpy > rp >  Vm > 2  cosð L Þ L > > > > > .. > < .

8 rp m  V 1  cosðrpL yÞ > L > > > rpy > rp >  Vm > 2  cosð L Þ L > > > > > .. > < .

9 > > > > > > > > > > > > > =

rpy > rp >  Vm > u  cosð L Þ > L > > > > > > > > > > . > > .. > > > > > > > > > > : rp m rpy ;  V  cosð Þ nþ1 L L

 rp L

Vu rpyu :  Vm  u cosð L Þ

ð60Þ

This relationship will be written in the following form

Vi ¼

rp Vm Vm i  cosð L  yi Þ i  V ¼  f1  V u ; u m V u  cosðrLp  yu Þ Vm u

ð61:aÞ

Uu

(a)

(b)

Fig. 11. Finite element mesh and constrained nodes: (a) for distortional and global buckling modes a longitudinal unknown degree of freedom is chosen and (b) for local buckling modes modes a transverse unknown degree of freedom is chosen.

1992

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

where yi is the longitudinal coordinate corresponding to node i. These are the constraint equations for the longitudinal displacements of the member. In a similar way the constraint equations for the transverse degrees of freedom are obtained L  U m  sinðrLp  yi Þ Um i  Vu ¼ m U i ¼ rp m i  f2  V u ; rp V u  cosð L  yu Þ Vu

ð61:bÞ

L  W m  sinðrLp  yi Þ Wm i  V ¼   f2  V u ; W i ¼ rp m i u V u  cosðrLp  yu Þ Vm u

ð61:cÞ

L  hm  sinðrLp  yi Þ hm  V u ¼ im  f2  V u : hi ¼ rp m i rp V u  cosð L  yu Þ Vu

ð61:dÞ

These equations are applied when the mesh is constrained to distortional or global GBT deformation modes. They cannot be used to constrain it to local modes, because in this case the longitudinal displacements are zero (Vu = 0). A transverse degree of freedom is chosen to be the unknown degree of freedom (see Fig. 11b), and the constraint equations change to

V i ¼ 0; Ui ¼

rp Um Um i  sinð L  yi Þ  U u ¼ im  f3  U u ; m rp Uu U u  sinð L  yu Þ

Wi ¼

rp Wm Wm i  sinð L  yi Þ  U u ¼ mi  f3  U u ; rp Um  sinð  y Þ Uu u u L

hm  sinðrLp  yi Þ hm  U u ¼ im  f3  U u : hi ¼ im rp U u  sinð L  yu Þ Uu

ð62:aÞ ð62:bÞ

ð62:cÞ

8 rp m1 9  V 1  cosðrpL yÞ > > L > > > > > > > rpy > < rp  V m1 = 2  cosð L Þ V u;1 L VðyÞ ¼ C m1   rp m1 . > >  V u;1  cosðrpLyu Þ > > L > > .. > > > > : ; rpy rp  V m1 nþ1  cosð L Þ L 8 rp m2 9  V 1  cosðrpL yÞ > > L > > > > > > > rpy > < rp  V m2 = 2  cosð L Þ V u;1 L þ C m2   rp m2 rpy . > > > > L  V u;1  cosð L u Þ > > .. > > > > : rp m2 ;  V nþ1  cosðrpL yÞ 8 L m1 9 8 m2 93 2 > > V1 > V1 > > > > > > > > >7 > > > 6 m2 > > > > > < V m1 < = =7 cosðrpyÞ 6C V ð1  C Þ 2 2 m1 7 6 m1 L   V u;1 ¼ 6 m1  þ 7 m2 . . > 6V u;1 > cosðrpLyu Þ V u;1 > .. > > .. > > >7 > > > > 5 4 > > > > > > > > : m1 ; : m2 ; V nþ1  V nþ1 ð63Þ where Cm1 and Cm2 are participation factors, sub-index u, 1 denotes first unknown degree of freedom, super-index m1 denotes mode 1, and m2 mode 2. The first step is to choose the second unknown degree of freedom, Vu,2. Fig. 12 shows the unknown degrees of freedom used in the present investigation. When choosing these degrees of freedom, it is important to take into account that their constraint factor should be different from zero. Subsequently, the participation factor Cm1 is expressed in terms of the two unknown degrees of freedom. This can be done from equation (63) in a similar way as in (58)

" ð62:dÞ

If a transverse unknown degree of freedom is chosen for distortional and global buckling, it will not be necessary to define two different sets of constraint equations.

4.2. Constraint equations for combined GBT deformation modes In some analyses it is necessary to combine two or more GBT deformation modes to obtain the required result. For example, the ‘‘pure’’ torsional-flexural elastic buckling load of a member subjected to compression is determined combining two modes; or the pure local, distortional and global elastic buckling moments of a channel section subjected to bending about the symmetry axis should be determined by means of a combination of two or more modes. Eqs. (51) and (52) in Section 3 showed the procedure to be followed to constrain a mesh to a specific combination of modes. However, as already pointed out, this procedure cannot be applied in commercial finite element software such as ANSYS. This problem is solved in the present investigation by means of constraint equations, in a similar way as in the previous section. The new constraint equations link all constrained node displacements to a set of unknown constrained degrees of freedom. There are as many unknown degrees of freedom as the number of combined modes required in the analysis. The derivation of these constraint equations is subsequently illustrated for the case of two mode combination. Only longitudinal degrees of freedom are included in the explanation. When two modes are combined, the longitudinal displacements can be expressed as a linear combination of two constraint equations such as (61.a)

V u;2 ¼

C m1

V m1 u;1



V m1 u;2

þ

ð1  C m1 Þ V m2 u;1

# 

V m2 u;2

 V u;1

ð64Þ

Constraint equations are determined by means of two unknown degrees of freedom that belong to the same cross-section. Thus, the cosinusoidal terms in (63) can be removed. Participation factor Cm1 is determined from (64) m2 V m2 u;2 =V u;1 m1 m1 m2 ðV u;2 =V u;1  V u;2 =V m2 u;1 Þ

C m1 ¼

¼ C1 þ C2 

þ

1 V u;2  m1 m2 m2 V u;1 ðV m1 =V  V =V Þ u;2 u;1 u;2 u;1

V u;2 V u;1

ð65Þ

Any longitudinal degree of freedom, Vi, located in yu can be expressed in a similar way as Vu,2 in (64)

"

Vi ¼

C m1

V m1 u;1



V m1 i

þ

ð1  C m1 Þ V m2 u;1

#



V m2 i

 V u;1

ð66Þ

Substituting (65) in (66), the constraint equations for the nodes of cross-section yu are obtained

Vi ¼

V m1 i V m1 u;1



V m2 i V m2 u;1

!

"

 C 2  V u;2 þ

V m2 i V m2 u;1

þ

V m1 i V m1 u;1



V m2 i V m2 u;1

!

#

 C 1  V u;1 ð67Þ

Eq. (67) can be extended to all longitudinal degrees of freedom of the member by including the cosinusoidal terms ( ! " ! # ) V m1 V m2 V m2 V m1 V m2 i i i i i  C  C  V  f1 ; Vi ¼   V þ þ  2 u;2 1 u;1 V m1 V m2 V m2 V m1 V m2 u;1 u;1 u;1 u;1 u;1

ð68Þ where f1 was defined in Eq. (61.a).The same expression of Cm1 can be used for the transverse degrees of freedom. As a consequence, their constraint equations can be directly set. For instance, the equation corresponding to the transverse rotation is

1993

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

( hi ¼

hm1 i V m1 u;1



hm2 i V m2 u;1

!

"  C 2  V u;2 þ

hm2 i

V m2 u;1

þ

hm1 i V m1 u;1



hm2 i V m2 u;1

!

#  C 1  V u;1

)  f2 ; ð69Þ

where f2 was defined in Eq. (61.b). Constraint equations for local buckling modes can be derived in a similar way. Vu2

5. Illustrative example Vu1

Fig. 12. Unknown constrained degrees of freedom for two mode combination.

E=210000 N/mm ν =0.3

2

E=200000 N/mm ν =0.3

1.67

20

13.1

89.7

3.33 100

2

68.8

125

(a) S1 (Schardt)

E=200000 N/mm ν =0.3

(b) S2 (Davies and Leach)

2

2.00

1.00

10

90

200

2

E=210000 N/mm ν =0.3

15

60

This section presents the calculation of the pure elastic buckling loads of a uniformly compressed channel member. The analysis is carried out on the cross-section shown in Fig. 13c, that was calculated in [21] by means of GBT based formulas. Loads are determined for the buckling of the member in one-half sine wave (r = 1 in (61)). Fig. 14a shows the finite element mesh and boundary conditions of a 600 mm long member. Only half of the member is included in the model. Symmetry boundary conditions are applied at one end, while at the other end a simple support is considered. The linear buckling analyses are performed in two steps. In the first step, a uniformly distributed normal pressure is applied at the hinged end of the member, and a linear analysis is carried out. This step is performed because ANSYS needs to calculate the internal stresses of the mesh to construct the geometric stiffness matrix of (50). In ANSYS, this preliminary calculation is not automatically carried out when performing a buckling analysis. In the second step, once the geometric matrix can be determined, ANSYS solves the eigenvalue problem. Between these two steps, constraint equations (61) (and (62)) are applied to the model to force the mesh to buckle in the chosen deformation mode. Table 1 contains the constraint factors used by the constraint equations for the member analysed in this example. The constrained model is shown in Fig. 14b, where it can be observed that there are 15 constrained nodes per constrained cross-section. Fig. 15 shows one of the pure buckling modes obtained (distortional), and Fig. 16 the elastic buckling loads for different halfwave lengths. In this figure, L1, L2 and L3 are the buckling loads of the first three local modes shown in Fig. 2, SD and AD are the symmetric and antisymmetric distortional buckling loads, and BX, BZ and T are the buckling loads of the global deformation modes of bending about X and Z axes, and torsion, respectively.

80

6. Verification of the constraining procedure

(c) S3 (Silvestre and Camotim)

(d) S4 (Ádány and Schafer)

Fig. 13. Geometry and material properties of the cross-sections analysed (dimensions in mm).

6.1. Individual GBT deformation modes Fig. 13 shows the cross-sections analysed to verify the accuracy of the constraining procedure. They are extracted from papers

Fig. 14. Finite element model: (a) mesh and boundary conditions. (b) GBT constraints.

1994

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

Table 1 Constraint displacement factors of cross-section S3.

U1 U2 U3 U4 U5 U6 U7 U8 U9 U10 U11 U12 U13 U14 U15 V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 V11 V12 V13 V14 V15 W1 W2 W3 W4 W5 W6 W7 W8 W9 W10 W11 W12 W13 W14 W15 h1 h2 h3 h4 h5 h6 h7 h8 h9 h10 h11 h12 h13 h14 h15

L1

L2

L3

SD

AD

BX

BZ

T

0.1208 0.0000 0.0000 0.0000 0.0000 0.0000 0.4777 0.7172 0.4777 0.0000 0.0000 0.0000 0.0000 0.0000 0.1208 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.1654 0.2512 0.2028 0.0000 0.0000 0.0000 0.0000 0.0000 0.2028 0.2512 0.1654 0.0000 0.0000 0.0121 0.0120 0.0091 0.0017 0.0083 0.0186 0.0193 0.0000 0.0193 0.0186 0.0083 0.0017 0.0091 0.0120 0.0121

0.5909 0.0000 0.0000 0.0000 0.0000 0.0000 0.3884 0.0000 0.3884 0.0000 0.0000 0.0000 0.0000 0.0000 0.5909 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.7628 1.0159 0.6436 0.0000 0.0000 0.0000 0.0000 0.0000 0.6436 1.0159 0.7628 0.0000 0.0000 0.0595 0.0583 0.0375 0.0053 0.0401 0.0376 0.0034 0.0242 0.0034 0.0376 0.0401 0.0053 0.0375 0.0583 0.0595

0.5543 0.0000 0.0000 0.0000 0.0000 0.0000 0.2149 0.5415 0.2149 0.0000 0.0000 0.0000 0.0000 0.0000 0.5543 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.6952 0.8642 0.4657 0.0000 0.0000 0.0000 0.0000 0.0000 0.4657 0.8642 0.6952 0.0000 0.0000 0.0559 0.0545 0.0322 0.0103 0.0367 0.0156 0.0219 0.0000 0.0219 0.0156 0.0367 0.0103 0.0322 0.0545 0.0559

3.8607 0.5767 0.5767 0.5767 0.5767 0.5767 3.9017 5.0101 3.9017 0.5767 0.5767 0.5767 0.5767 0.5767 3.8607 0.6804 0.1449 0.1031 0.0613 0.0195 0.0223 0.0223 0.0223 0.0223 0.0223 0.0195 0.0613 0.1031 0.1449 0.6804 17.0769 17.0769 12.1919 7.5533 3.4072 0.0000 0.0000 0.0000 0.0000 0.0000 3.4072 7.5533 12.1919 17.0769 17.0769 0.3284 0.3284 0.3202 0.2956 0.2545 0.1970 0.0985 0.0000 0.0985 0.1970 0.2545 0.2956 0.3202 0.3284 0.3284

4.7723 0.8864 0.8864 0.8864 0.8864 0.8864 1.5915 0.0000 1.5915 0.8864 0.8864 0.8864 0.8864 0.8864 4.7723 0.6660 0.1808 0.1166 0.0523 0.0119 0.0762 0.0381 0.0000 0.0381 0.0762 0.0119 0.0523 0.1166 0.1808 0.6660 17.5214 17.5214 11.7775 6.5440 2.3313 0.3504 0.3504 0.3504 0.3504 0.3504 2.3313 6.5440 11.7775 17.5214 17.5214 0.3886 0.3886 0.3716 0.3205 0.2355 0.1164 0.0367 0.0877 0.0367 0.1164 0.2355 0.3205 0.3716 0.3886 0.3886

0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.2273 0.2922 0.2922 0.2922 0.2922 0.2922 0.1461 0.0000 0.1461 0.2922 0.2922 0.2922 0.2922 0.2922 0.2273 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 1.3442 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000

2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 2.1005 0.3973 0.3973 0.2450 0.0927 0.0596 0.2119 0.2119 0.2119 0.2119 0.2119 0.0596 0.0927 0.2450 0.3973 0.3973 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000

1.5999 2.0573 2.0573 2.0573 2.0573 2.0573 1.0286 0.0000 1.0286 2.0573 2.0573 2.0573 2.0573 2.0573 1.5999 0.5063 0.3101 0.1609 0.0118 0.1373 0.2864 0.1432 0.0000 0.1432 0.2864 0.1373 0.0118 0.1609 0.3101 0.5063 4.0612 4.0612 3.3751 2.6890 2.0030 1.3171 1.3171 1.3171 1.3171 1.3171 2.0030 2.6890 3.3751 4.0612 4.0612 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457 0.0457

written by some of the main researches of GBT: cross-section S1 was analysed by Schardt in [5]; cross-section S2 by Davies and Leach in [22]; as mentioned above, cross-section S3 is extracted from a paper by Silvestre and Camotim [21]; and cross-section S4 was analysed by Ádány and Schafer in [11]. Nowadays, the MATLAB program developed in the present investigation to calculate the constraint factors does not allow for constraining all nodes of the cross-section. This is the reason why only the three configurations of cross-section constraints shown in Fig. 17 have been fully analysed herein. It will be shown,

however, that acceptable results are obtained with these configurations. The accuracy of the constraining procedure is verified comparing the results of the FEM analyses to the results of the GBTUL program. Only configurations b and c are considered for local buckling modes, since they are highly determined by the displacements of the intermediate nodes. The verification of the constraining procedure is carried out taking into account the three local modes shown in Fig. 2. Fig. 18 shows the ratios between the FEM local loads and the loads obtained by means GBTUL. It can be concluded that at

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

Fig. 15. Pure distortional buckling mode of cross-section S3.

1995

Fig. 19. Nb FEM vs Nb GBTUL for distorsional modes.

S3

Fig. 16. Curves of pure buckling loads of cross-section S3. Fig. 20. Nb FEM vs Nb GBTUL for global modes.

Fig. 17. Number of constrained nodes per cross-section.

Fig. 18. Nb FEM vs Nb GBTUL for local modes.

least fifteen constrained nodes per cross-section are needed to get accurate results. Distortional modes (modes SD and AD in Fig. 2) are first calculated constraining only 6 nodes per cross-section (Fig. 17a). The re-

sults are good for members longer than the critical distortional length, that is, the length corresponding to the smallest distortional buckling load (LcrD). On the contrary, the distortional buckling loads of short members are rather underestimated. Results improve as the number of constraints per cross-section is increased. This is verified by repeating the analyses with 9 and 15 constrained nodes. It is observed, however, that more nodes should be constrained to obtain accurate values for the shortest members (L < 0.3  LcrD).fig. 19. For global buckling modes, nine constrained nodes per crosssection are enough to get acceptable results (Fig. 20). A small sensitivity analysis of the results to the number of constrained cross-sections is also carried out. S1 and S4 are analysed with different number of cross-sections linked to the unknown displacement. Fig. 21 shows the result of the analysis. It can be observed that the number of constrained cross-sections can be reduced to 10 per half-sine wave with a small degradation of results. It is important to underline that reducing the number of constraint equations will reduce the calculation time in ANSYS. On the other hand, it can also be observed that highly constrained meshes lead to overestimated buckling loads of distortional and global modes. These values are obtained when the ratio between the spacing of constrained sections and the element size equals one (spacing/element size = 1), that is, when all crosssections of the mesh are constrained. A similar problem is found when working with cFSM. In this method, constraints provoke an increase of rigidity and, consequently, the distortional and global critical loads obtained are slightly higher than the GBT loads. A complete discussion on this issue can be found in [11–13]. cFSM and GBT results coincide only when the Poisson’s ratio is zero.

1996

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

Fig. 22. Results of the constrained analyses: (a) S1, (b) S2, and (c) S4.

Fig. 21. FEM/GBTUL ratio vs number of constrained sections per half- wave: (a) local, (b) distortional and (c) global buckling loads.

The same occurs in the constrained models used in the present investigation. If all the cross-sections of the mesh are constrained, the accuracy of FEM results is good only when the Poisson ratio is zero (or very small) in both FEM and GBT models. Constrained finite element analyses were performed using the same number of constraints but different spacing/element size ratio. It was observed that when this ratio equals 2 the loads obtained are rather accurate. Fig. 22 shows the curves of elastic buckling loads of cross-sections S1, S2 and S4. They are calculated constraining 15 nodes per cross-section for local and distortional modes, and 9 nodes per cross-section for global modes. In the longitudinal direction, more than 10 cross-sections per half-sine wave are constrained, and a spacing/element size ratio equal to 2 is applied. In Figs. 23 and 24, the results of the constraining procedure are compared

to the results of GBTUL and CUFSM. It can be observed that the FEM buckling loads are rather similar to GBTUL loads, but for the distortional mode of very short members. When compared to CUFSM values, the differences for the loads of the distortional mode are more significant. This is due to the mentioned effects of cFSM on the rigidity of the member. Similar differences should be observed for the global buckling loads. However, Fig. 24c does not show them because, in the end, FEM global loads were compared to values obtained by means of the classical global buckling analytical formulas.

6.2. Combined GBT deformation modes Fig. 25 shows some results of constrained analyses carried out combining GBT modes. A compressed S3 cross-section is considered in the calculations, and the chosen combinations are: L1 + SD, L2 + AD, T + BX and T + BX + AD. The results of the analyses are compared to GBT buckling loads of the same mode combinations in Fig. 26. A rather good accuracy is achieved.

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

1997

Fig. 23. Nb FEM vs Nb GBTUL: (a) local, (b) distortional, and (c) global. Fig. 24. Nb FEM vs Nb CUFSM: (a) local, (b) distortional, and (c) global (analytical).

On the other hand, Figs. 27–29 show the results of the constrained analyses carried out on cross-section S3 subjected to bending about the x axis (BX). The combined modes are: L1 + L2 + L3 for pure local mode; SD + AD for pure distortional mode (see Fig. 27); and T + BZ for pure global mode. Generally speaking, the buckling moments obtained are acceptable except for short members buckling in the pure local mode. The FEM local buckling moments are higher than the GBT local buckling moments because the number of combined modes is too small. If the FEM values are compared to GBT moments determined with GBTUL considering the same number of deformation modes, the results become acceptable (see Fig. 29b). Therefore, it can be concluded that in this case more than three GBT deformation modes should be combined to obtain accurate results for pure local buckling moments. Calculations with combined modes were carried out on the other cross-sections of Fig. 13. It was verified that the accuracy achieved in these analyses was similar to the one shown by S3. In view of the results shown in this section, it seems that the proposed procedure for the calculation of pure buckling loads works rather well when applied to members subjected to pure compression or bending. Further investigations should be per-

Fig. 25. Nb for combined modes of S3 cross-section.

formed to validate it for members with different loading conditions, such as variable normal stresses along the longitudinal axis or beams under different distributed loads.

1998

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

Fig. 26. Nb FEM vs Nb GBTUL for combined modes of S3.

Fig. 27. Pure distortional buckling mode of cross-section S3 subjected to bending. Fig. 29. Mb FEM vs Mb GBTUL for S3: (a) FEM local loads compared to GBTUL local loads calculated combining all local modes; (b) FEM local loads compared to GBTUL local loads calculated combining only the first three local modes.

Fig. 28. Pure buckling moments of cross-section S3. Fig. 30. Local (L1) and symmetric distortional buckling loads of C–C S3.

7. Members with general boundary conditions All calculations shown in previous sections were performed considering simple supports at both ends of the member (S–S). In the present section, calculations applying two additional boundary conditions are included. The boundary conditions considered are: C–C, both ends clamped; and C–S, clamped and simply supported. The aim of this section is to show that the proposed procedure can easily be applied to members with general boundary conditions. The sinusoidal function used for w(y) in equation (4) has to be changed to take into account the effect of the two new end conditions [23,24]

C—C : wðyÞ ¼ sin

C—S : wðyÞ ¼ sin

r  p  y  p  y  sin ; L L 

 r  p  y ðr þ 1Þ  p  y þ ðr þ 1=rÞ  sin : L L

ð70Þ

ð71Þ

As a consequence, parameters f1, f2 and f3 in constraint equations (61) and (62) also change (calculations are made for pure modes, no combinations). The new parameters can be calculated from (70) or (71), depending on the boundary condition applied

1999

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

Fig. 31. Local (L2) and antisymmetric distortional buckling loads of C–C S3.

Fig. 34. Nb FEM vs Nb GBTUL of S3 with C–S boundary conditions.

clamped–clamped S3 member are drawn with dashed lines in Figs. 30 and 31. Calculations are performed for one to four half-waves (r = 1, 2, 3, and 4). It should be noticed that the critical mode for a specific length of a member must be expressed by means of Eqs. (70) and (71) combining several half-wave terms (r = 1, 2, 3, ...) [23]. For instance

wðyÞ ¼

nhw X

sin

 r  p  y

r¼1

L

 sin

p  y : L

ð73Þ

In the present investigation, however, this approach is not applied to obtain the critical modes. Instead, parameters f1, f2 and f3 are removed from constraint equations. The member is left free to buckle in the ‘‘longitudinal’’ direction, without any half-wave constraint. The idea is that the longitudinal shape of the buckling mode is a result of the linear buckling analysis, and it is just determined by the actual end boundary conditions introduced in the model. (This may have some advantages when boundary conditions are very complex, see [25]). New analyses are carried out removing the half-wave constraints. Each cross-section is independently constrained as shown in Fig. 32. The constraint equations become

Fig. 32. Independently constrained cross-sections.

Vi ¼

Vm i  V us ; Vm u

ð74:aÞ

Ui ¼

Um i  U us ; Um u

ð74:bÞ

Fig. 33. Nb FEM vs Nb GBTUL of S3 with C–C boundary conditions.

Wi ¼

0

f1 ¼

w ðyi Þ ; w0 ðyu Þ

ð72:aÞ hi ¼

wðy Þ f2 ¼ 0 i ; w ðyu Þ

ð72:bÞ

wðyi Þ : wðyu Þ

ð72:cÞ

f3 ¼

The four cross-sections shown in Fig. 13 are analysed with the new end supports. Finite element models are generated with C–C and C–S boundary conditions at member ends, and constrained according to the constraint equations that can be derived from parameters (71). The Nb that results from the analysis of a

Wm i  U us ; Um u

hm i  U us ; Um u

ð74:cÞ

ð74:dÞ

where Vus and Uus are the unknown degrees of freedom of cross-section s for the longitudinal and transverse degrees of freedom, respectively. In Eqs. (74.b–74.d), the transverse displacements are not expressed in terms of the longitudinal unknown displacement. The removal of the half-wave constraints provokes that no relationship can be set between these two types of degrees of freedom. Constraints (74) are applied to the four cross-sections analysed in the present investigation. Figs. 30 and 31 show the Nb curve obtained in the analysis of cross-section S3 (continuous line), and Figs. 33 and 34 the accuracy of results for the C–C and C–S bound-

2000

M. Casafont et al. / Computers and Structures 89 (2011) 1982–2000

ary conditions. Similar accuracy is obtained in the analyses of the other cross-sections investigated. 8. Conclusions One of the main features of GBT is that allows for calculating pure elastic buckling loads and modes. The aim of the present investigation has been to incorporate this feature into the Finite Element Method. The constrained Finite Strip Method has provided the key to success, since it already solved the problem. Therefore, it can be said that the main contribution of the paper is the translation of the GBT and cFSM concepts into the FEM. It is underlined that this contribution is just a tentative first step towards a constrained Finite Element Method. There are still many issues to be solved, as pointed out at the end of this section. The original formulation of the cFSM has been changed. It has been adapted for its use in commercial finite element software, such as ANSYS. Displacements are not directly constrained by means of transformation matrices. Constraint equations are used instead, which are, just from the formal point of view, of slightly more general application. Constraint equations are defined for different types of analyses: compressed members, members subjected to bending, members with different end supports, individual mode calculation or combined mode calculation. Therefore, it is demonstrated that the constraining procedure presented allows for the same calculations as the existing GBT and cFSM methods. The performance of the constraining procedure is evaluated by means of analyses on four different channel sections. Generally speaking, the results obtained are rather accurate. They are similar to the results of the GBTUL and CUFSM programs. Nowadays, the procedure is being applied with a limited number of constrained nodes per cross-section. In the near future, it will be extended to constrain all nodes of the cross-section. Then, better accuracy will be obtained, especially for distortional modes of short members. It will also be possible to constrain local buckling modes of higher order. Furthermore, constraining all nodes of the model will solve one of the problems of the procedure being applied: visual identification of the buckling mode is still needed. Non-constrained nodes provoke the presence of ‘‘spurious’’ modes among the results of the analysis. This is the reason why the modes that result from the analyses should be visually classified. On the other hand, it is also demonstrated that if all cross-sections of the mesh are constrained, slightly overestimated values of distortional and global buckling loads are obtained. Therefore, it seems that constraining all nodes will not solve all the problems, and that more investigation should be performed on this issue. Another point which needs further research is the role of the non constrained degrees of freedom hx and hz. It is necessary to study how to constrain them and verify the effect of these additional constraints on the final results of the linear buckling analyses. Finally, it is pointed out that the procedure presented herein will hardly replace GBT and cFSM in general applications. The cost of the combined GBT–FEM procedure in terms of model generation and, especially, computing time is high compared to the other methods. The method may have some advantages in specific applications. For instance, when analysing members with complex boundary conditions. However, as it is based on the beam theory assumptions of GBT, all the advantages of the shell finite element

models cannot be incorporated in the procedure. It is not possible to calculate pure buckling loads of perforated members, members with thickness variation along the length or members with any other discrete variation. References [1] Silvestre N, Camotim D. GBT buckling analysis of pultruded FRP lipped channel members. Comput Struct 2003;81:1889–904. [2] Basaglia C, Camotim D, Silvestre N. Global buckling analysis of plane and space thin-walled frames in the context of GBT. Thin-Walled Struct 2008;46:79–101. [3] Basaglia C, Camotim D, Silvestre N. GBT-based local, distortional and global analysis of thin-walled steel frames. Thin-Walled Struct 2009;47:1246–64. [4] Camotim D, Basaglia C, Silvestre N. GBT buckling analysis of thin walled steel frames: a state-of-the-art report. Thin-Walled Struct 2010;48:726–43. [5] Schardt R. Lateral torsional and distortional buckling of channel- and hatsections. J Constr Steel Res 1994;31:243–65. [6] Beregszászi Z, Ádány S. The effect of rounded corners of cold-formed steel members in the buckling analysis via the direct strength method. In: Topping BHV, Costa LF, Barros RC, editors. Proceedings of the 12th international conference on civil, structural and environmental engineering computing, CC2009. Funchal: Civil-Comp Press; 2009. [7] Joó AL, Ádány S. FEM-based approach for the stability design of thin-walled members using cFSM base functions. Period Polytech 2009;53:61–74. [8] Ádány S, Joó AL, Schafer BW. Buckling modes identification of thin-walled members using cFSM base functions. Thin-Walled Struct 2010;48:806–17. [9] Bebiano R, Pina P, Silvestre N, Camotim D. GBTUL – buckling and vibration analysis of thin-walled members. DECivil/IST, Technical University of Lisbon; 2008. [10] Ádány S, Schafer BW. Buckling mode decomposition of single-branched open cross-section members via finite strip method: derivation. Thin-Walled Struct 2006;44:563–84. [11] Ádány S, Schafer BW. Buckling mode decomposition of single-branched open cross-section members via finite strip method: application and examples. Thin-Walled Struct 2006;44:585–600. [12] Ádány S, Silvestre N, Schafer BW, Camotim D. Buckling analysis of unbranched thin-walled members using cFSM and GBT: a comparative study. In: Camotim D, Silvestre N, Dinis PB, editors. Proceedings of international colloquium on stability and ductility of steel structures, SDSS 2006. Lisbon: IST Press; 2006. [13] Ádány S, Silvestre N, Schafer BW, Camotim D. GBT and cFSM: two modal approaches to the buckling analysis of unbranched thin-walled members. Adv Steel Constr 2009;5:195–223. [14] Schafer BW, Ádány S. Buckling analysis of cold-formed steel members using CUFSM: conventional and constrained finite strip methods. In: Eighteenth international specialty conference on cold-formed steel structures, Orlando; 2006. [15] Casafont M, Marimon F, Pastor MM. Calculation of pure distortional elastic buckling loads of members subjected to compression via the finite element method. Thin-Walled Struct 2009;47:701–29. [16] Schardt R. Verallgemeinerte technische biegetheorie. Heidelberg: SpringerVerlag; 1989. [17] Silvestre N, Camotim D. First-order generalised beam theory for arbitrary orthotropic materials. Thin-Walled Struct 2002;40:755–89. [18] Ádány S. Buckling mode classification of members with open thin-walled cross-sections by using the Finite Strip Method. Research Report, Johns Hopkins University; 2004. [19] Silvestre N, Camotim D. Second-order generalised beam theory for arbitrary orthotropic materials. Thin-Walled Struct 2002;40:791–820. [20] Cook RD, Malkus DS, Plesha ME. Concepts and applications of finite element analysis. New York: John Wiley & Sons; 1989. [21] Silvestre N, Camotim D. Distortional buckling formulae for cold-formed steel C and Z-section members. Part II – Validation and application. Thin-Walled Struct 2004;42:1599–629. [22] Davies JM, Leach P. Some applications of generalized beam theory. In: Eleventh international specialty conference on cold-formed steel structures, St. Louis; 1992. [23] Silvestre N, Camotim D. Distortional buckling formulae for cold-formed steel C and Z-section members. Part I – Derivation. Thin-Walled Struct 2004;42:1567–97. [24] Li Z, Schafer BW. Finite strip stability solutions for general boundary conditions and the extension of the constrained finite strip method. In: Topping BHV, Costa LF, Barros RC, editors. Trends in civil and structural engineering computing. Stirlingshire: Saxe-Coburg Publications; 2009. [25] Camotim D, Silvestre N, Basaglia C, Bebiano R. GBT-based buckling analysis of thin-walled members with non-standard support conditions. Thin-Walled Struct 2008;46:800–15.