Constrained shell finite element method for elastic buckling analysis of thin-walled members

Thin-Walled Structures 145 (2019) 106409

Contents lists available at ScienceDirect

Thin-Walled Structures journal homepage: www.elsevier.com/locate/tws

Full length article

Constrained shell finite element method for elastic buckling analysis of thinwalled members

T

Sheng Jina,b,∗, Zhanjie Lic, Fang Huanga, Dan Gana,b, Rui Chenga,b, Gaofeng Dengd a

School of Civil Engineering, Chongqing Univ., Chongqing, 400045, China Key Laboratory of New Technology for Construction of Cities in Mountain Area (Chongqing University), Ministry of Education, Chongqing, 400045, China c Department of Engineering, SUNY Polytechnic Institute, Utica, 13502, USA d General Architectural Planning & Design Research Institute of Chongqing University Co., Ltd., China b

A R T I C LE I N FO

A B S T R A C T

Keywords: Constrained finite element method Constrained finite strip method Elastic buckling analysis Thin-walled member Distortional buckling

The objective of this paper is to develop a constrained shell Finite Element Method (cFEM) based on a force approach for elastic buckling analysis of thin-walled members. The new cFEM is able to separate the general deformation of thin-walled members into the three fundamental deformation mode classes, namely Global (G), Distortional (D), and Local (L), to enable the modal decomposition and identification. In this paper, four forcebased mechanical criteria are defined to separate these mode classes. These mechanical criteria are implemented with the general shell finite element formulation without any special treatment of the element formulation. The constraint matrices of the G, D, and L mode classes are then constructed. Numerical examples are presented to demonstrate the capabilities of the new cFEM in modal decomposition and identification. In particular, the modal decomposition and identification results are compared with the cFSM solutions. Applicability of the new cFEM to other shell FE formulation and different loading conditions are illustrated. All these numerical examples demonstrate the potential of the developed cFEM in taking advantages of the modeling capability of the existing shell FE method.

1. Introduction The elastic buckling analysis has played a critical role in understanding the behavior of a thin-walled member due to the high slenderness ratio of each individual element comprising of a thin-walled member. Additionally, the critical buckling stress usually serves as a good indicator of the ultimate strength of the thin-walled member. However, the buckling analysis of the thin-walled member can be quite challenging due to the many possible deformation modes that are generally encountered: local (plate), distortional, and global. As reflected in the current design specifications, e.g. AISI [1], each buckling mode is treated uniquely due to its significantly different post-buckling behavior. Furthermore, possible mode interaction among them can complicate the analysis even more. Therefore, the buckling analysis necessitates not only the traditional prediction of the critical stress but also an ability to identify buckling modes appropriately. Many numerical methods are available to perform the elastic buckling analysis of thin-walled members, such as Finite Element Method (FEM), Finite Strip Method (FSM) [2], and Generalized Beam Theory (GBT) [3–6]. Even though all of them are capable of predicting critical stress, not all

∗

of them are capable of mode identification and separation (i.e., FSM and FEM). This gap has been partially filled with the development of constrained Finite Strip Method (cFSM) [7] and more recently constrained Finite Element Method (cFEM) [8,9]. Including the GBT, all these methods can decompose the general buckling modes into the fundamental buckling mode classes, Global (G), Distortional (D), Local (L) and Shear and Transverse Extension (ST). A brief review of the development of these methods is provided in the following to gain a concept about the buckling modal analysis of thin-walled members. GBT is implemented in the form of an enriched beam element transformed from nodal to modal degrees of freedom (DOF) through the inclusion of the cross sectional deformations. Conventional GBT [3–6] did not consider the shear deformation and transverse extension, and the orthogonality among the basic modes are with respect to the warping stiffness C and transverse stiffness B [10]. Camotim and Silvestre have greatly expanded the GBT capabilities with a posteriori inclusion of particular deformation modes including the shear deformation modes [10–12] to a variety of applications on buckling analysis, such as open and closed thin-walled sections [13,14], members with orthotropic material [15,16], thin-walled frames [17,18], and

Corresponding author. School of Civil Engineering, Chongqing Univ., Chongqing, 400045, China. E-mail address: [email protected] (S. Jin).

https://doi.org/10.1016/j.tws.2019.106409 Received 19 April 2019; Received in revised form 20 August 2019; Accepted 16 September 2019 0263-8231/ © 2019 Elsevier Ltd. All rights reserved.

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

Δm, Pm

Δm = Ux ∪ Us, Pm=Nx ∪ Ns, the membrane subset of the universal set of ΔE and PE (loaded DOFs in the Global and Distortional deformations, also constrained DOFs in the Local deformation, see Fig. 2) Ux, Nx translations and forces of all the nodes in the X-direction (longitudinal direction) Us, Ns translations and forces of all the nodes in their transverse tangential directions Δb, Pb Δb = Δ¯m , Pb = P¯m , the local-plate bending subset of the universal set of ΔE and PE (unloaded DOFs in the Global and Distortional deformations, also unconstrained DOFs in the Local deformation)

Notation n k p nm nb ΔE

PE

number of cross-sections (including the end sections) number of nodes on a cross-section number of intersection nodes on a cross-section number of membrane DOFs number of bending DOFs displacement vector, including three translation DOFs (U) each and three rotation DOFs (θ) each of all the nodes (see Fig. 1) load vector corresponding to the ΔE, including three forces (N) and three moments (M) on each node

implementation based on the general shell FE formulation is desired. It was shown in Ref. [34] that the characteristics of the applied loads can be used to classify the deformation modes of thin-walled members and limited applications in modal identification were shown. Hence, inspired by the works of Adany [8,9], Khezri and Rasmussen [31,32], the objective of this paper is to develop an alternative cFEM with expanded mechanical criteria of defining the fundamental deformation modes based on the general shell FE formulation (e.g., shell elements in ANSYS). The new mechanical criteria, significantly different from the mechanical definitions in current GBT [5,15], cFSM [20], and cFEM [8,9], utilize the stiffness matrix of the first-order elastic analysis solution to specify the deformation modes. The implementation also requires no distinction between main and subordinate nodes, which paves the way to the mode class definitions that could be consistent for both the prismatic and curved cross-sections. In addition, the new cFEM takes advantage of the general shell FE formulation and does not require a special treatment of element formulation. The force-based criteria enable the orthogonality between the basic modes with respect to the stiffness. Furthermore, the method can take advantage of the existing computational codes such as ANSYS to perform the buckling analysis and hence potentially nonlinear analysis in the future. The paper begins with the shell finite element representation by explaining the degrees of freedom (DOFs), the node definitions, and the decomposition of the displacement field and the stiffness matrix. The next section of the paper focuses on the deformation mode definition – the four mechanical criteria used for defining each mode class. Then the derivation of constraint matrices for modal decomposition and identification are presented in the following section. Finally, a variety of numerical examples related to modal decomposition and identification is performed to highlight the feasibility of the method. Applicability to other shell FE formulation and loading conditions are demonstrated. The paper finishes with a discussion of limitations, and potential improvements and future research possibilities.

composite thin-walled members [15]. cFSM leverages the efficiency of FSM [19] by incorporating the same GBT mechanical assumptions of the buckling modes. By constraining the strip's membrane and plate bending behaviors, constraint matrices are constructed to separate the degrees of freedom into those consistent with L, D, G, and other (e.g., Shear and Transverse extension - ST) deformation modes [20]. Most recently, Li et al. [19]. provides a review of the development and applications of cFSM including the extension of the constraint techniques for spline FSM [21,22]. Meanwhile, given the limitations of cFSM applicable to prismatic members with constant moment gradient and specific end boundary conditions, significant extensions are undertaken in the domain of the shell FEM by applying the cFSM base functions or GBT cross-sectional deformation modes [23,24]. One direction in these efforts is to utilize the cFSM base function to approximately identify the mode participations of the buckling modes [25] or even the general deformation from the collapse analysis of thin-walled members [26]. Another effort is to apply the constraint of the GBT deformation modes into the shell FEM models in order to calculate the pure buckling loads from shell FEM [27,28]. More recently, Ádány [8,9] developed a new cFEM based on a novelly formulated shell FE and then applying the same mechanical assumptions as those in cFSM. On the other hand, from the formulation point of view, constructing the constraint matrices varies significantly depending on the particular features of the methods although they all share the same underlying fundamental assumptions for the buckling mode definitions (note, true for cFSM and GBT as well). New implementation techniques have been proposed, which significantly differ from the traditional GBT implementation in Ref. [29] or cFSM in Ref. [30] for other application purposes or for better efficiency. Notably, Khezri and Rasmussen [31,32] proposed an energy-based formulation for determining pure modes of buckling for prismatic thin-walled members. Piccardo et al. [33] utilizing the Kantorovich's semi-variational method developed a quadratic function based on the energies for transversal bending, longitudinal bending and axial strains for formulating the fundamental buckling modes that can be used in the GBT analysis. Jin [34] introduced a new force-based approach to separate the fundamental buckling modes based on the characteristics of the applied and resultant loads. In addition, Becque proposed new modal decomposition concepts based on equivalent nodal forces [35] or the extremal participation of normalized bending and membrane strain energies [36]. In general, compared to others, the development of the cFEM bears numerous merits in the potential applications by exploiting the versatility of the modeling capability of the shell finite element method to handle complex geometry, loading, and boundary conditions. This provides the opportunity for modal decomposition and identification similar to cFSM that can replace a laborious and subjective procedure employing visual investigation to classify the buckling modes in elastic buckling analysis. The only realization of the cFEM is developed by Ádány [8,9] for buckling analysis of thin-walled members, however, the shell FE formulation is specially defined. A more generalized

2. Shell finite element representation The method proposed herein utilizes the existing shell FE formulations such as those in the commercial package ANSYS. The core of the shell FE formulation is intact and only the manipulation of the integrated stiffness matrix is performed as explained in this section. A shell FE model of a thin-walled member modelled in ANSYS is shown in Fig. 3 as an illustration. The global coordinate system is denoted as X–Y–Z, with the X-axis parallel to the longitudinal direction of the member. For the cross section, a local Cartesian coordinate system x-sw applies, with x- and s-axes parallel to X-axis and along the mid-line, respectively, as depicted in Fig. 4. Fig. 4 also illustrates the specification of the nodes: intersection node and ordinary node. No necessary distinction further than that on the main nodes as in current cFSM is needed. Fig. 1 illustrates the six DOFs of a node in a shell finite element. The membrane stresses of the element are mainly related to the longitudinal 2

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

PE* = K e⋅ΔE*

(3)

By rearranging the components in the matrix and vectors according to the classification of the membrane DOFs and the bending DOFs, Eq. (3) is rewritten as * T * ⎧ Pm ⎫ ⎡ Kmm Kbm ⎤⋅⎧ Δm ⎫ = ⎥ ⎨ ⎨ Pb* ⎬ ⎢ K Kbb ⎦ ⎩ Δb* ⎬ ⎭ ⎩ ⎭ ⎣ bm

(4)

3. Deformation mode definition The fundamental deformation mode classes: Global (G), Distortional (D), and Local (L), are defined through four mechanical criteria explained in this section. It is worth noting that in current cFSM and cFEM, the other mode classes, Shear (S) and Transverse Extension (T), are usually defined as well as separate mode classes. In the proposed method, they are embedded in the G and D mode classes.

Fig. 1. The six DOFs of a node.

translation Ux and the transverse tangential translation Us. Therefore, for each ordinary node (nodes other than the intersection node), there is only one transverse tangential translation Us, while for each intersection node, both the two transverse translation DOFs, Uy and Uz, should be considered because both will produce membrane strains in at least one adjacent wall as shown in Fig. 2. Note any transverse translation of an intersection node could also create some bending in at least one adjacent wall, while according to GBT, it is exactly this sort of bending that should be excluded from the Local mode bending and forms the Distortional bending. As a result, for a member with n cross-sections considered along the member length, and with k nodes (including p intersection nodes) for each cross section, the total number of the membrane stress correlating DOFs is nm = n(2k + p). In this paper hereinafter, all these DOFs are termed as membrane DOFs and denoted as Δm. It is worth noting that theoretically the independent rotational DOF about the normal to the shell surfaces, i.e. the drilling DOF θw in Fig. 1, should be considered as one of the membrane DOFs, though they have no stiffness associated based on the shape functions of most shell elements [37]. A relatively small stiffness is added in the shell FEM formulation based on different measures (see e.g. Ref. [37]) to prevent numerical instability or to achieve better performance. The authors found that the inclusion of the drilling DOF into the membrane DOFs has a negligible effect on the results but rise more complication of the constraint matrix's derivation. Hence, in this study, the drilling DOF is excluded from the membrane DOFs Δm. Meanwhile, all the DOFs other than membrane DOFs are termed as bending DOFs (Δb) in this paper. The total number of the bending DOFs is nb = n(4k-p), including the three rotational DOFs of all the nk nodes and the one normal translation each for all the n(k-p) ordinary nodes. For a thin-walled member modelled in shell FEM, the load-displacement equation of a linear elastic analysis problem is usually given by

PE0* = K e0⋅ΔE0*

3.1. Definition of local modes Criterion #1: For L deformation modes, all the displacements of membrane DOFs are null. This criterion is consistent with GBT, in which local modes are characterized by the absence of in-plane displacements of the intersection nodes and no involvement of the cross-sectional longitudinal displacements [15]. Hence, Criterion #1 leads to the following expression:

ΔmL = 0

(5)

3.2. Properties of the applied loads on the global and distortional modes Criterion #2: For G and D deformation modes, the applied loads corresponding to the bending DOFs are null. This means that

PbG = 0, and PbD = 0

(6)

Thus the orthogonality of the local deformation to the other deformations can be obtained with: T

T

T

{ΔLΕ} ·PEG(D) + {ΔLm} ·PmG(D) + {ΔLb } ·PbG(D) = 0

However, this criterion is not enough to separate the G/D modes. From the force equilibrium point of view, the transverse deformation of D modes according to GBT can be regarded as the result of transverse forces acting on the intersection nodes. All these forces are acting on the membrane DOFs. Meanwhile, the G modes, generally known as bending,

(1)

which is defined according to the global coordinate system X-Y-Z. Where ΔE0 is a column vector consists of the six DOFs (Ux, Uy, Uz, θx, θy and θz) of all the nodes, PE0 is the corresponding load vector, and Ke0 is the stiffness matrix. In this paper, the Global, Distortional, and Local deformation modes are expressed by the superscripts ‘G’, ‘D’ and ‘L’ separately. The superscript ‘*’ in Eq. (1) stands for a general category of G, D, and L under consideration. The transverse translations Uy, Uz, and the corresponding loads, Ny and Nz of all the ordinary nodes, are transformed into translations and forces according to the local coordinates, by T N U ⎛ Ns ⎞ = ⎛ cos α − sin α ⎞ ⋅⎛ y ⎞ , and ⎛ y ⎞ = ⎛ cos α − sin α ⎞⋅⎛ Us ⎞ ⎝ Nw ⎠ ⎝ sin α cos α ⎠ ⎝ Nz ⎠ ⎝Uz ⎠ ⎝ sin α cos α ⎠ ⎝Uw ⎠ ⎜

⎟

⎜

(7)

⎟

(2) where α is the inclined angle from the Y-axis to the local s-axis on the positive X-plane. Thus Eq. (1) is transformed into

Fig. 2. The membrane DOFs on a cross section. 3

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

represents the entrywise product between two matrices of the same dimensions. edt and edw are arbitrary 3 × 1 and 4 × 1 vectors, respectively. 3.4. Additional loading property of the distortional modes The orthogonality between the G/D modes,

{PED }T⋅ΔEG

=

∑ [{ eNxD }T⋅e UxG + { eNsD }T⋅e UsG] = 0 e

(10)

leads to the conditions the applied loads on each cross section ‘e’ should satisfy:

{ e NsD }T⋅[cos α s, sin α s, ys ∘ cos α s-z s ∘ cos α s ] = {0 0 0}

(11)

and

{

e

NxD }T⋅[ 1k × 1, e

y, z , ω ] = {0 0 0 0}

(12)

e

where, Ns and Nx are the column vectors consisted of all the transverse tangential forces and all the longitudinal forces applied on the nodes of the cross section ‘e’, respectively. Though obtained mathematically, Eqs. (11) and (12) do hold specific mechanical meanings, which can be summarised as: Criterion #4: For D deformation modes, all the resultant forces, moments, and bimoment of the applied loads acting on each cross section are null. Although the Criterion #4 is defined through properties of the external applied loads, it does not deviate from the traditional concepts of the distortional mode. In GBT, the orthogonality condition between G/ D modes writes [38]:

Fig. 3. FEM model of a deformed lipped channel.

∫ω u xD⋅u xG⋅dΩ = 0

(13)

where Ω is the area of the cross section. Based on simplified constitutive relations, it is deduced in Ref. [5] for the distortional mode that i) the cross-sectional tangential stresses are in self-equilibrium, and ii) the profile of the cross-sectional normal membrane stresses is in accordance with that of the primary warping, which means Fig. 4. Local coordinate and Node specifications.

∫ω σxD⋅u xG⋅dΩ = 0

(14)

torsion, and axial deformation of a member, can be correlated directly to the seven resultant forces, moments, and bimoment acting on each cross section along the member length. Based on these, additional criteria related to G/D modes can be derived in the following subsections.

Thus in the original concepts of the distortional mode, the resultant internal forces, moments and bimoments on a cross section are all null.

3.3. Additional displacement properties of the global deformation modes

The proposed criteria allow the possible decomposition of the general deformation into the fundamental G, D, and L mode classes. This leads to the constraint matrices associated with these fundamental mode classes similar to the cFSM. Thus, the modal decomposition (i.e., calculating the critical loads for pure buckling modes) and modal identification (i.e., identifying the mode participations of a general deformation mode into G, D, and L) are enabled for shell FEM. In this section, the constraint matrix for each deformation mode class will be constructed by applying the above criteria through establishing the interrelationships of the considered mode class within all the DOFs. Then these constraint matrices can be introduced into the general eigenvalue equations for modal decomposition and identification.

4. Constraint matrix for shell finite element

Criterion #3: For global deformation modes, i) the “plane crosssection” assumption and/or torsional warping mode are applicable to the longitudinal membrane DOFs, and ii) the “rigid-body” assumption is applicable to the transverse membrane DOFs. The assumptions of plane cross-section/torsional warping and rigidcontour are the traditional hypotheses for the bending and torsional deformation modes. Since Criterion #2 has defined a condition for each bending DOF, this criterion excludes conditions related to these DOFs. Based on this criterion, the associated DOFs can be written as the following for each cross section ‘e’ along the member length as: e

UsG = [cos α s, sin α s, ys ∘ sin α s −z s ∘ cos α s ]⋅e d t

e

UxG

(8)

4.1. DOF interrelationship of local deformation modes

and

= [ 1k × 1, y, z , ω ]⋅e d w

The DOF vector can be written as: (9)

L

Δ ΔEL = ⎧ mL ⎫ = RL ⋅ΔbL ⎨ Δb ⎬ ⎩ ⎭

where αs is the inclined angle from the Y- axis to each Us direction of a cross section, ys and zs are the Y- and Z-coordinates of the node corresponding to each Us DOF, respectively, and y, z, and ω are respectively the Y-, Z-, and sectorial coordinates of each node. The sign ‘∘’

L

(15)

where R is the constraint matrix of the local deformation modes. Because the local deformation modes satisfy Eq. (5) – Criterion #1, that 4

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

J ⋅ΔmD = 0

leads to

0nm × nb ⎤ RL = ⎡ ⎢ In b ⎥ ⎣ ⎦

where, (16)

T −1 J = T T⋅(Kmm − Kbm ⋅Kbb ⋅Kbm )

According to Criterion #2, Eq. (6) leads to

= Kbm·ΔmG +Kbb⋅ΔbG = 0 From Eq. (17),

ΔbG

=

ΔbG

(17)

can be calculated as

−1 −Kbb Kbm⋅ΔmG

(18)

ΔmG

consisting of the longitudinal translation Ux and the transverse with tangential translation Us, the displacement property equations (8) and (9) of all cross sections of Criterion #3 can form the ΔmG as the following:

ΔmG

1 G 1 ⎧ Us ⎫ ⎧ dt ⎫ ⎪ ⋮ ⎪ ⎪ ⋮ ⎪ ⎪ n G⎪ ⎪ G ⎪ ⎪ nd t ⎪ ⎧Us ⎫ ⎪ Us ⎪ = T ⋅d = ⋅ = T = ⎨ 1d w ⎬ ⎨UxG ⎬ ⎨ 1UxG ⎬ ⎩ ⎭ ⎪ ⎪ ⋮ ⎪ ⎪ ⎪n ⎪ ⎪ ⋮ ⎪ ⎪ dw⎪ ⎪ nUxG ⎪ ⎭ ⎩ ⎭ ⎩

D

⎧ Δm1 ⎫ = 0 [ J1 J2 ]⋅ D ⎬ ⎨ Δm2 ⎩ ⎭

D D Δm1 =J3⋅Δm2

(19)

Inm ⎤⋅T =⎡ −1 ⎢− Kbb ·Kbm ⎥ ⎣ ⎦

Should be mentioned that sometimes the J1 matrix might be nearly singular, in this case MATLAB can still perform this calculation, while the result may be inaccurate — elements with extremely high absolute values can be found in the resulted J3 — these shall introduce remarkable errors into the following processes. An automatic optimization measure can solve this problem: Considering that the i j-th element possesses the largest absolute value in J3, and which is extremely high, this means that the i-th DOF in D D Δm1 is not a good choice, and the j-th DOF of Δm2 might be better. We can exchange them and recalculate the J3 through Eq. (30) based on the updated J1 and J2, or more efficiently, to reassign the J3 elements by

(20)

[J3]kl =

According to Criterion #2, Eq. (6) leads to the following equation similar to Eqs. (17) and (18) for the G deformation modes:

=

−1 −Kbb ·Kbm⋅ΔmD

(22)

=0

k = i, l ≠ j k = i, l = j

(31)

(32)

D

where, R is the constraint matrix of the distortional deformation modes as the following,

Inm ⎤·X ⋅⎡ J3 ⎤ RD=⎡ D −1 ⎢− Kbb ⎢ ·Kbm ⎥ ⎦ ⎣ Inm − 7n + 6 ⎥ ⎣ ⎦

23)

Meanwhile, according to Criterion #4, combining equations (11) and (12) of all the cross sections together, leads to:

T T⋅PmD

⎨− [J3]i l /[J3]i j ⎪ J ⎩1/[ 3]i j

D ΔED = RD⋅Δm2

Then applying Eq. (4) to the D mode class leads to: −1 T PmD = (Kmm − Kbm ⋅Kbb ⋅Kbm )⋅ΔmD

⎧[J3]k l − [J3]k j ⋅[J3]i l /[J3]i j k ≠ i , l ≠ j ⎪[J3]k j /[J3]i j k ≠ i, l = j

The formula (31) is a matrix multiplication operation to the J3 and computationally much efficient. Repeating the DOF exchanging process D D between the Δm1 and Δm2 enough times will lead to a reasonable choice D of the Δm1 DOFs and the resulted J3. Finally, combining Eqs. (22) and (29) into the overall DOFs for D D : modes, ΔED , it can be written based on Δm2

(21)

4.3. DOF interrelationship of distortional deformation modes

ΔbD

(30)

J3 = −J1 \J2

where RG is the constraint matrix of the global deformation modes as the following,

RG

(29)

where the (7n-6) × (nm-7n+6) matrix J3 can be obtained using the matrix backslash operator in MATLAB:

G

⎧ Δm ⎫ = R G ⋅d ⎨ ΔbG ⎬ ⎩ ⎭

(28)

where the 7n × (7n-6) matrix J1 and the 7n × (nm -7n+6) matrix J2 are D D and Δm2 , respectively. the J's columns corresponding to Δm1 D Thus the Δm1 can be calculated, through:

where d is an arbitrary 7n-element column vector, and the matrix T is nm × 7n. Eqs. (18) and (19) indicate that the G modes are defined with 7 DOFs for each cross section. This is different from the conventional 4 in GBT as the differential relationship between the longitudinal and transverse displacements is not integrated here, thus the corresponding shear deformation is not separated from the G modes as GBT and cFSM do. According to Eqs. (18) and (19), the DOF interrelationships of the G modes are written as:

ΔEG =

(27)

Note that although the matrix J has 7n rows, the rank of J is actually (7n-6). Because the PmD obtained by Eq. (23) is self-equilibrated, which means in all the 7n null resultant force conditions defined in Eq. (26), six equilibrium conditions are satisfied automatically if the other (7n-6) conditions have been satisfied. Knowing that the Eq. (26) defines (7n-6) independent conditions for the ΔmD, (7n-6) elements in the ΔmD can be determined by the other elements: arbitrarily select (7n-6) elements in the ΔmD , designate them D D and the others in ΔmD as Δm2 , Eq. (26) can be rewritten as as Δm1

4.2. DOF interrelationship of global deformation modes

PbG

(26)

(33)

and the nm × nm matrix XD represents the row exchanging operation D D and Δm2 from ΔmD , i.e. used to separate the Δm1

(24)

D

Δ ⎫ ΔmD =XD⋅⎧ m1 D ⎬ ⎨ Δm2 ⎩ ⎭

where, D

N PmD = ⎧ sD ⎫ ⎨ Nx ⎬ ⎩ ⎭

(25)

(34)

4.4. Buckling mode decomposition

Eq. (24) defines seven conditions (of null resultant forces) each for all the n cross sections. Then, introducing Eq. (23) into Eq. (24), leads to

Need to mention that the basic modes discussed above are indeed 5

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

and constrain them as 0s. Then, a constraint equation consisting of enough boundary conditions thus eliminating the redundant rigid-body deformations can be written:

modes of deformation and independent of buckling problems. While these deformation modes serve mainly for the distinguishing between different buckling types, like in the GBT and cFSM. Buckling deformations and critical loads of a thin-walled member are obtained by solving the general eigenvalue problem:

K e · ϕ = Λ · K g ⋅ϕ

D ΔED = RD′⋅Δm2_II

where the Δm2_II is the redundant components of Δm2 by removing the six components of Δm2_I, and RD’ is the resultant matrix by removing the six corresponding columns. Now Eq. (39) can be rewritten in a more specific format with a fully spanned basis as

(35)

where, Ke is the stiffness matrix, Kg is the geometric stiffness matrix, vectors of ϕ represent the buckled shapes, and Λ is a diagonal matrix of eigenvalues. It is assumed in buckling mode decomposition problem that the member buckles according to a specific pure deformation mode. For modal decomposition, the general eigenvalue problem is constrained to the mode class through the constraint matrices constructed above, thus in a similar fashion like cFSM, this leads to

K e∗⋅ϕ∗

=

Λ∗ ⋅Kg*⋅ϕ∗

L

⎧c ⎫ ΔEO = [ RL RG RD′]⋅ c G ⎨ D′⎬ ⎩c ⎭

ΔE* = R*⋅c *

where, the superscript ‘*’ may be replaced by ‘ ’, ‘ ’, or ‘ ’ to express the decomposed eigenvalue problem of G, D or L buckling modes (or any combination thereof), Λ* is a diagonal matrix of eigenvalues, ϕ* is the buckling modes in column, and the elastic and geometric stiffness matrices of the constrained buckling problems are

K e* = [R*]T ·K e⋅R*

D

L

(42)

It is worth noting that different choices of the Δm2_I DOFs, i.e. different systems of statically determinate constraints introduced, may lead to different D sub-deformations ΔED from the same original deformation ΔEO , yet these D sub-deformations only differ from each other for a rigid-body shift and have the same value of strain energy. The strain energy corresponding to each mode class as well as the total strain energy of the original deformation can be written as:

(37)

and

Kg* = [R*]T ·K g⋅R*

(41)

According to the solutions of Eq. (41), the separated L, G and D subdeformations of the original deformation ΔEOcan be written:

(36) G

(40)

E∗ =

(38)

Should be noted that the “pure” buckling mode might be different from any real buckling (i.e. general buckling) mode of the member. That is because the buckling deformation is somewhat assumed and, accordingly, the equilibrium is not between the hypothesis loads * * PEg (=λ*·Kg⋅ΔE* ) and the elastic forces PEe (=K e⋅ΔE* ) but their equivalent * and the [R*]T ·PEe *. loads, i.e. the [R*]T ·PEg

1 *T (ΔE ) ⋅K ⋅ΔE∗ 2

(43)

and

EO

= EL + E G + ED

(44)

where ‘*’ represents the mode class, G, D, and L. The mode participation factor of each mode class, according to the portion of strain energy, is * cEnergy = E ∗/ E O

(45)

4.5. Buckling mode identification 5. Numerical examples A general deformation mode can be categorized into the fundamental deformation mode class, G, D, and L in terms of their participations. This is called modal identification. The constraint matrices of the mode classes defined above are able to span the original shell FEM's DOF space and provide an alternate basis organized in G, D, and L subspaces. In Ref. [34], the modal identification was explored through the strain energy contributions though the constraint matrices were defined differently. A modal identification process similar to cFSM is presented here. However, it should be noted that the current constraint matrices do not have a separate Shear and Transverse extension (ST) subspace. The ST modes in this paper are included inside the GD modes. The constraint matrices RL, RG, and RD defined above forms the basis R needed to span the original shell FEM DOF space. For an arbitrary deformation ΔEO , their participations can be obtained through a linear combination of the column vectors of the basis R, i.e. solving the c in the equation

ΔEO = R⋅c

To demonstrate the capability of the developed cFEM, a set of numerical examples on modal decomposition and modal identification are given. The results are specifically compared and discussed with the currently available solutions using cFSM. The studied section is a Steel Stud Manufacturers Association (SSMA) stud section SSMA 800S250–97 [39] under simply supported axial compression. This section was the same one used for the illustrative example of cFSM in Ref. [20]. A material with the elastic module of 203.5 GPa and Poisson ratio of 0.3 is assumed and the round corners of the cross section are neglected. The shell FEM solution is conducted in ANSYS using the 4nodes, 24-DOF SHELL63 (KEYOPT(3) = 2 [37]) element. The cross section is discretized by 3, 7, and 1 intermediate nodes in the flanges, web, and lips, respectively. The longitudinal element dimension is about 20 mm. The boundary conditions are set to simulate the simply supported boundary conditions in cFSM: a static determinate boundary condition (i.e. constraining only 6 nodal DOFs of the whole shell FE modelled member) is applied when calculating the pre-buckling stresses, while constraining the Y- and Z-translations of all the nodes on both ends and the X-translation of an arbitrarily selected node when performing the linear buckling analysis, as explained in Ref. [34]. An alternative method for the settings of the boundary conditions can be seen in Ref. [40].

(39)

where c is a column vector that includes the contribution coefficients for each vector inside each deformation class. However, for the constraint matrices RL, RG, and RD defined above, the total column vectors are nm+nb+6, which exceeds the total DOFs of the shell FEM, nm + nb, by 6. This is because the rigid-body displacements (6 DOFs) are accommodated both in global and distortional modes defined above. Any arbitrarily defined statically determinate constraints can be specified to remove this redundancy [34]. This redundancy removal can be performed on either RG or RD. Taking the RD for example, for the constraint Eq. (32), pick out six DOFs, Δm2_I, that can form a statically determinate constraint system,

5.1. Modal decomposition After establishing the shell FE model, the elastic stiffness matrix is extracted from ANSYS, and a MATLAB program is utilized to construct the constraint matrices (16), (21) and (33) based on the elastic stiffness 6

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

and developed cFEM is a little larger compared to the differences between L and G modes. The D mode defined herein does not separate the ST modes, thus lower pure D curves compared to cFSM shall be obtained in most cases, like the one shown in Fig. 15, while in this example a higher D curve than cFSM is obtained. Moreover, when compared to the D+ST mode of cFSM, the difference is around 5% with the cFEM solution higher. This is mainly due to another difference between the two mechanical definitions of the D modes: the null requirement of the applied normal forces is enforced for the thin-plates here, while enforced only for the ‘equivalent multi-span beam’ of the cross-section in cFSM. This will result in differences due to the discrepancies between the bending behaviors of beams and plates. The authors found that for the buckling cases when the coupling between the L and D modes are high like the current example (indicated by no distortional minimum in the FSM signature curve), the discrepancy caused by the second effect of the D definitions might lead to higher D curves compared to the cFSM as shown here. A comparison to the GBT results has also been conducted. The GBT analytical results of the L, D, D+ST, G, and G+Shear are essentially in agreement with the cFSM results, thus similar conclusions when compared to this work can be drawn.

matrix. Then, Eqs. (15), (20) and (32) are respectively introduced into the shell FE model in the form of ANSYS constraint equations. Finally, performing the ANSYS linear buckling analysis under these constraints obtains the solutions of Eq. (36), i.e. the decomposed buckling results. The decomposition is performed for each mode class, G, D, and L as shown in Fig. 5, in which the critical stress of 1st buckling mode is plotted as a function of the physical length of the column. The three pure mode solutions of the same section were simulated using the cFSM for simply supported boundary conditions under axial loading as well using CUFSM's general boundary conditions with 1–40 longitudinal terms [41]. In addition, the signature curve of this section under axial loading (i.e., only m = 1 longitudinal term) is also plotted in Fig. 5 and the horizontal axis should represent the half-wave length instead. Note, the cross sectional level discretization is the same between the FEM and CUFSM models. As shown in Fig. 5, the pure mode solutions of G, D, and L from cFSM and the developed cFEM agree well in general. The selected buckling modes in the figure that clearly demonstrate the corresponding mode classes of G, D, and L also highlight the correctness of the mode definition mechanism proposed herein. The longitudinal waves of the buckling modes with the changing member length are also accurately captured. In particular, the pure Local buckling deformation agrees extremely well between the two methods. This is attributed to the exact same definition of the local modes in the two methods. Even though slight difference does exist but that is due to the difference between the shell finite element and finite strip formulations – the FSM element is in general slightly stiffer than the FEM element [42]. As is well known that the shear strains do affect the flexural and torsional buckling behaviors, while these are not considered in the conventional G mode. As can be seen from Fig. 5, when shear deformation is considered, the G+Shear buckling results of cFSM is slightly smaller than its pure G results. The corresponding shear considered in this work leads a G curve agree well with the cFSM G+Shear curve . Meanwhile, the difference between the developed cFEM and regular FSM (roughly 10% higher) shares the similar mechanical discrepancy as discussed in Ref. [7] for cFSM and regular FSM due to the constraint of transverse extension. The difference of the pure Distortional solutions between the cFSM

5.2. Modal identification The general buckling analysis solutions can be extracted from ANSYS. With the constraint matrices established in the same way as those for modal decomposition in the previous section, the modal identification can be performed by calculating the participations in terms of their strain energy contributions following the proposed approach in the previous section. The same computational model used for modal decomposition in the last sub-section is adopted here to illustrate the modal identification capability of the developed cFEM. 5.2.1. Modal participation of the 1st buckling mode The lowest critical stresses of the studied member are plotted in Fig. 6a as a function of the member length against the regular FSM solutions for general S-S boundary conditions with 1–40 longitudinal terms (numbers of half-waves). At the same time, the mode

Fig. 5. Buckling mode decomposition of a lipped channel under axial compression. 7

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

Fig. 6. Modal identification of the member (the 1st buckling mode).

5.2.2. Modal participation of higher-order buckling modes To further highlight the modal identification capability of the developed cFEM, the examination of higher modes at a fixed length is provided in this sub-section. This is more aligned with the engineer's design convention: model the physical member and identify the associated buckling modes of G, D, and L in higher modes for the critical stresses. To illustrate this, the first 30 buckling modes of the member with a length of 1400 mm are extracted from the ANSYS's eigen-buckling solutions. The identified G, D, and L participations and selected buckling shapes are summarised in Table 1 and Fig. 7, respectively. In addition, a comparison with the FSM results of this member is attempted here. The FSM solution was calculated using CUFSM, the same as those in the previous section. The reason this is only an attempt is because the eigen solution varies with different eigen solvers. Even though the 1st eigen solution may be close or the same, the higher modes do not necessarily align with each other order by order, which makes the direct comparison hard. Visual inspection has been applied

participations of G, D, and L as functions of the member length are calculated using the proposed modal identification approach and plotted in Fig. 6b against the cFSM participations as well. It is worth noting that the modal identification of cFSM utilizes the orthogonal basis with a strain energy norm and the summation of the participation factors adopts the L-2 norm (see Ref. [43] for details). Fig. 6b shows that the modal identification results of the two methods consistently agree well with each other over the studied length with an only slight difference. As can be seen, when the member length is smaller than 2500 mm, the lowest critical loads are all corresponding to a local mode with a dominant local participation (i.e., > 80%). Beyond that length, the member is dominated by the global mode (i.e., Euler buckling mode). It is worth noting that unlike cFSM, the transverse extension modes were not separated in the developed method here, which shall be an immediate research task. Though not the focus here, Fig. 6a does highlight the excellent agreement of the stability solutions of the FSM and FEM.

Table 1 GDL participations for the first 30 buckling modes. this paper # Mode 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

CUFSM

σcr (MPa)

G%

D%

L%

168.4 170.3 170.7 176.4 177.2 184.9 189.7 195.7 205.3 208.6 216.7 223.3 225.9 239.8 257.8

0.2 0.3 0.2 0.1 0.4 0.1 0.6 0.1 0.9 0.1 1.3 0.1 2.1 0.1 0.0

2.1 3.6 1.3 0.8 6.7 0.5 13.2 0.3 26.5 0.2 48.5 0.2 69.6 0.1 0.1

97.7 96.1 98.6 99.1 92.9 99.4 86.3 99.6 72.6 99.7 50.2 99.8 28.3 99.8 99.9

Δ#

+1 −1

+1 −1

this paper

Δσcr (MPa)

Δ G%

Δ D%

Δ L%

+2.3 +2.1 +2.5 +2.8 +1.8 +3.1 +1.6 +3.4 +1.3 +3.8 +1.1 +4.2 +0.8 +4.7 +5.2

+1.0 +1.1 +0.9 +0.8 +1.3 +0.7 +1.5 +0.6 +1.6 +0.6 +1.5 +0.5 +1.3 +0.5 +0.5

+0.5 +0.8 +0.3 +0.2 +1.3 +0.1 +2.4 +0.1 +4.1 +0.1 +4.9 +0.0 +3.2 +0.0 +0.0

−2.9 −3.8 −2.4 −2.0 −5.0 −1.7 −7.1 −1.5 −9.8 −1.3 −11.1 −1.2 −8.4 −1.1 −1.0

# Mode 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

8

CUFSM

σcr (MPa)

G%

D%

L%

277.3 285.5 298.3 320.7 344.4 369.4 387.4 395.8 423.4 452.3 480.3 482.4 513.8 546.0 546.4

0.0 5.8 0.0 0.0 0.0 0.0 66.4 0.0 0.0 0.0 92.1 0.0 0.0 0.1 0.0

0.1 78.1 0.1 0.1 0.0 0.0 28.5 0.0 0.0 0.0 7.8 0.0 0.0 0.5 0.0

99.9 16.1 99.9 99.9 99.9 99.9 5.1 99.9 100.0 100.0 0.1 100.0 100.0 99.4 100.0

Δ#

Δσcr (MPa)

Δ G%

Δ D%

Δ L%

+2 −1

+5.8 +0.5 +6.5 +7.2 +8.0 +8.9 +0.4 +9.8 +10.9 +11.9 +1.1 +13.1 +14.3 +19.2 +15.6

+0.4 +1.8 +0.4 +0.4 +0.4 +0.4 +1.6 +0.4 +0.4 +0.4 −1.4 +0.4 +0.4 +0.9 +0.4

+0.0 −0.1 +0.0 −0.0 −0.0 −0.0 −9.5 −0.0 −0.0 −0.0 −7.5 −0.0 −0.0 +0.2 −0.0

−1.0 −4.3 −0.9 −0.9 −0.8 −0.8 −0.4 −0.8 −0.8 −0.8 +0.2 −0.8 −0.7 −2.3 −0.7

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

Fig. 7. Comparison of the modal identification for high order buckling deformations.

modes are coupled by their participations. For instance, if that dot has color more close to Red meaning that the buckling mode has large local mode participation. From Fig. 8, it is interesting to notice that all these higher modes can potentially form similar signature curves with higher modes and different longitudinal terms from CUFSM as shown in Ref. [19]. Fig. 8 also shows the transverse symmetric, asymmetric, and tri-symmetric patterns of the buckling modes among the higher modes as indicated in the figure. Meanwhile, with the help of the rendered colors, it is easy to identify the variation of the modal participations (indicating the dominant mode: RGB → LGD) with respect to the member length. If all the L dominated (i.e. the L participation is the largest among the three), D dominated, and G dominated dots in Fig. 8 are chosen, respectively, and their lower boundary lines are drawn, the so-called separated buckling curves in this paper are obtained and shown in Fig. 9. These curves should represent the associated buckling modes of L, D, and G as an engineer would identify for the design purpose by examining the higher modes. However, the difference is that in Fig. 9,

to align the modes and Table 1 shows the differences of the critical stress and participations compared to CUFSM solutions for all the 30 modes. Even though the differences are small for most of them, the attention should be probably given to those more closely matched, for instances, those shown in Fig. 7. These mode shapes are visually identical and their critical stresses and participations are all in good agreement with respect the two methods. 5.2.3. Another perspective view of buckling mode participations Fig. 6a only shows the lowest critical stress curve as a function of the member length from stability solutions. At any particular member length, the higher modes can be examined and mode participations can be identified as demonstrated in the last sub-section. If the member length is limited from 50 mm to 3200 mm and the critical stress capped at 1525 MPa, for each member length, the higher modes fall into the cap can be put into Fig. 8 as dots and theirs mode participations are L indicated by the color with a RGB format (Red: Green: Blue = cEnergy : D G cEnergy : cEnergy ). Hence, the mixture of the color indicates how the LGD

Fig. 8. Color rendered modal participations of the FE linear buckling modes. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.) 9

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

Fig. 9. Separated buckling curves based on modal identification.

All the previous examples are performed with the SHELL63 element. For the same member and same mesh, the pure mode solutions and mode participations are recalculated using a 4-node, 24-DOF shell element SHELL181 element in ANSYS. SHELL181 is more prevalent used because of its plastic capability and suitability for analyzing thin to moderately-thick shell structures [37]. The decomposed pure mode solutions using different shell FE are shown in Fig. 10. Slight difference of the critical stresses exists, especially for local buckling mode, resulting from the different formulation of the element types; however, the difference is very consistent along the length, indicating that the good applicability of the method to a general shell FE type. Meanwhile, the separated buckling curves of G, D, and L using SHELL181 FE similar to the one in the last section using SHELL63 FE are shown in Fig. 11. Though the critical stresses slightly differ but the consistent agreement of the curves indicates again the developed cFEM can adopt any general shell FE - no special formulation is required.

the identification is done through the rigorous modal identification with quantification merit other than the subjective visual inspection. In Fig. 9, the pure G, D, and L solutions from the developed cFEM and signature curve from the FSM are also shown for comparison. According to AISI [1], the first minimum and the final descending branch of the FSM signature curve are identified, respectively, as the local and global buckling stresses and utilized in the Direct Strength Method (DSM) [44] for the estimation of the loading capacity. Can be seen that the proposed separated L and G curves provide the same critical results. The distortional buckling in this example is difficult to identify through the FSM signature curve due to the indistinct second minimum. While the separated D curve provides multiple minimums with the same value of 217.4Mpa, which is quit close to the analytical solution of the distortional buckling, σcrd = 232.4Mpa, calculated using the (Eq. 2.3.1.3-2) of AISI S100-16w [1]. One word should be mentioned is that the minimum value of the separated D relies on the judgement rule of D dominant. In this paper it is simply prescribed that a higher D participation than the L and G participations means a D dominant buckling, while this needs a further discussion but out the scope of this work. The regularities of the separated D and L curves indicate that for some members, only one linear buckling analysis may provide reliable and adequate critical stresses required by DSM. Take the 1400 mm length member for example, the one general FEM buckling analysis results in σcrl = 168.4 MPa(#1 mode, with 97.7% L participation), σcrd = 225.9 MPa(#13 mode, with 69.6% D participation), and σcre = 387.4 MPa (#22 mode, with 66.4% G participation), as listed in Table 1. These SHELL FE numerical elastic buckling solutions, automatically identified but justified by the visual inspects in the last subsection, should be permitted to be utilized in DSM.

7. Applicability to closed and branched sections A branched cross-section with a closed cell as depicted in Fig. 12a is considered in this example. The application of a generalized cFSM to this XHL cross-section can be seen in Ref. [45], In the current paper the

6. Sensitivity to shell finite element type One of the merits of the developed cFEM is the utilization of the generally formulated shell FE. More specifically, no specially formulated shell FE is needed and no decoupling of the membrane and bending behaviors is necessary. Moreover, the developed cFEM can directly adopt any existing formulation of the shell FE. Hence, to demonstrate the applicability of the developed cFEM and highlight the sensitivities of the element type, the modal decomposition solutions are illustrated through the same member used above with a different element type.

Fig. 10. Modal decomposition – different element type. 10

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

loading, boundary condition, and material are set the same as in Ref. [45], i.e. the member is axially compressed, simple-simple supported, and E = 210 GPa. The only exception is that a more normally used Poisson's ratio of 0.3 but not 0 is assumed here. As for the FE model, SHELL63 is used, 32 nodes and 32 elements are considered for the cross-section, as depicted in Fig. 12b, and the longitudinal element dimension is set to be about 6 mm. Though the cross-section contains branches and a closed cell, the definitions of the basic modes remain as before, and the modal decomposition and identification can be implemented without any change. The identifications on the general buckling of shell FE modelled thin-walled member with different lengths are depicted in Fig. 13. As before, the Red, Green, and Blue color participations represent the L, G and D modal participations, and help recognizing the “curves” corresponding to different longitudinal half-wave numbers and different transverse deform patterns. Considering the lower boundary of the “curves” with single longitudinal half-wave, they agree well with the FSM signature plotted in Fig. 14. The first three buckling shapes of the 25 mm length member are plotted in Fig. 13. In these three deformations, the transverse translations of all the 8 intersection nodes, and thus the transverse tangential translations of all the 16 thin-walls, are negligible, thus they are identified as L dominant. With the development of the half-wave length, the transverse translations of the intersection nodes become remarkable, thus the 1st mode changes to torsional (G) one, while the 2nd and 3rd modes become Distortional ones (with the critical half-wave lengths of about 175 mm), and the 3rd mode comes to flexural (G) at last. Buckling mode decomposition results are depicted in Fig. 14. The pure G curve and buckling shapes are in good accordance with the FEM 1st linear buckling results for members longer than 1000 mm, and the pure D curve and the buckling shape well agree with the lowest D dominated general buckling results. The pure L buckling shape is more like the 2nd but not the 1st L dominated general buckling plotted in Fig. 13, a little bit deviate from the expectation but does not seriously affect the resulted critical stress and critical half-wave length of the L

Fig. 11. Separated buckling curves – different element type.

Fig. 12. An XHL cross-section's (a) Geometry and (b) Discretization.

Fig. 13. Buckling mode identification of an XHL column. 11

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

Fig. 14. Buckling mode decomposition of an XHL column.

Second, as expected, the moment gradient provides the benefits of the buckling strength for global buckling as well as distortional buckling where the member length is relatively long. In addition, the modal decomposition of D modes of non-uniform bending also captures the descending trend with the beam length, which is different from the uniform bending case. The relative localized D mode (also L mode when the beam is relatively long) for non-uniform bending beam is also expected as shown in Fig. 16. Third, when the beam is short, the critical stresses of pure L modes under non-uniform bending are significantly lower than those of the uniform bending. This is due to the shear effect for short beams due to the moment gradient. The developed cFEM by utilizing the analysis capacity in ANSYS that can handle shear correctly captures this effect. The buckling mode is a local mode due to the shear stress as shown in Fig. 16. The shear effect diminishes when the beam length increases. The moment gradient gradually shows some beneficial effect for local buckling (i.e., higher critical stress) though small. The mode participations of the buckling modes with the lowest four critical stresses are analyzed using the developed modal identification approach. The rendered RGB colored dots corresponding to all these linear buckling modes are plotted in Fig. 17 against the pure L, D, and G

mode. This is because the torsion of the closed-cell is strictly forbidden by Criterion #1 for the pure L mode. 8. Applicability to other loading cases One of the merits of the developed method is to utilize the modeling capacity of the shell FEM to handle complex boundary conditions and loading cases. To illustrate this, the same member used before are modified with loading cases to highlight the feasibility of the developed cFEM. In this example, two loading cases of a simply supported beam are considered to highlight the capability of the developed cFEM. The first beam loading case is uniform bending by applying concentrated end moments about the cross-sectional major axis (i.e., uniform moment along the member). The other one is non-uniform bending. The concentrated end moment is applied on one end only, therefore the moment varies along the member length. One issue should be mentioned is about the settings of the boundary conditions. When performing the linear buckling analyses, as before, transverse translations of all the nodes on both ending sections and an arbitrarily selected longitudinal translation are constrained. When calculating the pre-buckling stresses, for the uniform bending case, a static determinate BC (6 DOF constraints) is applied which, like in the axially compressed case, will results in pretty “pure” 1st stresses as predefined in GBT and FSM. However, when non-uniform bending is applied, the forces produced in these 6 constraints of the static determinate BC should cause serious stress concentrations in the member. Thus for the non-uniform bending case, the applied BC when calculating the prebuckling stresses is the same as that used in buckling analysis. For comparison, the pure mode solutions of cFSM and general FSM buckling results, all with general Simple-Simple boundary conditions, of the uniform bending beam are also calculated using CUFSM. Note CUFSM doesn't have the capability to consider the loading cases with moment gradient. The modal decomposition results of the uniform bending beam are plotted in Figs. 15 and 16 against the CUFSM solutions and the non-uniform bending beam result, respectively. First, for uniform bending, the CUFSM solutions showed three buckling regimes: local, distortional, and global as the beam length varies from short to long. The critical stresses of pure L, D, and G modes from the developed cFEM agree well with CUFSM solutions in their respective regime.

Fig. 15. Buckling mode decomposition of the beam subjected to uniform bending. 12

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

questionable. Meanwhile, the L modes deforming between the main nodes will be significantly stiffened for fine discretization needed for the curved geometry. If the definition of the main and subordinate nodes is discarded, the force-based approach herein holds advantages to quickly implement the needed criteria based on the similar definitions of the membrane and bending DOFs in separating the modes. The development of cFEM for curved sections is underway. Ideally, this development may also address the current issue of prismatic members with rounded corners. Finally, in this paper, the new cFEM is implemented in the commercial software package ANSYS. All the analysis capabilities of ANSYS can be potentially pursued including the constrained post-buckling and nonlinear collapse analyses. More research is needed to expand this capability.

10. Conclusions Fig. 16. Buckling mode decomposition of beams subjected to uniform and nonuniform bending.

The development of a new constrained shell Finite Element Method (cFEM) paves a way for a broader application of the shell finite element method by enabling the modal decomposition and identification in elastic buckling analysis of thin-walled members. In this paper, four force-based mechanical criteria are proposed to separate the Global, Distortional, and Local modes. The new cFEM takes advantage of the general shell FE formulation and does not require a special treatment of element formulation. The force-based criteria enable the orthogonality between basic modes with respect to stiffness. In addition, the implementation also requires no distinction between main and subordinate nodes, which paves the way to the mode class definitions that could be consistent for both the prismatic and curved cross-sections. The construction of the constraint matrices for the mode classes follows the proposed four force-based criteria. The numerical examples of modal decomposition and identification of an axial compression member demonstrate excellent agreement with cFSM results though small differences do exist. It is noted that the differences stem from not only the difference between FEM and FSM models but also the slight discrepancy in the definitions of the GD modes between the new cFEM and cFSM. Furthermore, the modal decomposition example using another shell FE formulation highlights one of the merits of the new cFEM: the independence of the shell FE formulation. Meanwhile, numerical examples on other loading cases demonstrate another merit of the method in the potential of exploiting the superior modeling capability of shell FEM. In particular, the new cFEM is able to capture the shear effect for short beams and correctly identify and decompose the shear buckling modes. Finally, the limitations and future developments to more general thin-walled cross section are also identified. The potential of the new cFEM for post-buckling and nonlinear collapse analyses can be foreseen but need more research.

mode buckling curves. The variation of the color from red to blue to green matches the trends of the pure modes well and indicates the local, distortional, and global regimes as well as the potential transition regions where the modes are greatly coupled.

9. Discussions The mechanical definitions of the buckling modes in the developed cFEM herein do not address the separation of the shear and transverse extension modes (ST) that are commonly defined in GBT and cFSM. These ST modes are currently embedded in the GD modes. For most cases, the contribution of ST modes to the modal identification is small, thus they will have negligible influence on the results except for certain nonlinear analyses as demonstrated in GBT [46] and cFSM [26] (e.g., slightly larger than 2%). Nonetheless, the separation of the transverse extension mode from GD is an immediate task of the authors. In addition, it should be mentioned that the difference demonstrated in this study between the cFSM and cFEM integrates the difference of mesh patterns, stiffness functions, and types and numbers of DOF considered as a whole. The implementation difference of the forcebased criteria in this paper compared to the existing cFSM cannot be fully identified unless the developed definitions are implemented on the same FSM model. Future research is warranted in this aspect. Moreover, flat elements are often used to discretize thin-walled members with curved cross-sections, such as the circular cylindrical tube. The challenge exists how to define the nodes used to discretize the curved section. If defined as main nodes, with no transverse translation permitted at any main node for the L modes the D modes will overwhelmingly dominate the deformation space, which is very

Fig. 17. Buckling mode identification vs. pure GDL buckling for non-uniform bending. 13

Thin-Walled Structures 145 (2019) 106409

S. Jin, et al.

Acknowledgments [23]

The support of this work by the National Program on Key Research and Development Project (No. 2017YFC0703805), the Natural Science Foundation Project of CQ CSTC (No. cstc2017jcyjAX0341), 111 Project of China (No. B18062) and the Fundamental Research Funds for the Central Universities (No. 106112015CDJXY200001) is gratefully acknowledged.

[26]

References

[27]

[24] [25]

[1] American Iron and Steel Institute AISI S100-16, North American Specification for the Design of Cold-Formed Steel Structural Members, 2016 Edition, American Iron and Steel Institute, Washington (DC), 2016. [2] Y.K. Cheung, The finite strip method in the analys of elastic plates with two opposite simply supported ends, Proc. Inst. Civ. Eng. 40 (1) (1968) 1–7. [3] R. Schardt, Verallgemeinerte Technische Biegetheorie [Generalised Beam Theory], Springer Verlag, Berlin, 1989 ([in German]). [4] R. Schardt, Generalized beam theory-an adequate method for coupled stability problems, Thin-Walled Struct. 19 (2–4) (1994) 161–180. [5] J.M. Davies, P. Leach, First-order generalised beam theory, J. Constr. Steel Res. 31 (2–3) (1994) 187–220. [6] J.M. Davies, P. Leach, D. Heinz, Second-order generalised beam theory, J. Constr. Steel Res. 31 (2–3) (1994) 221–241. [7] S. Ádány, Buckling Mode Classification of Members with Open Thin-Walled CrossSections by Using the Finite Strip Method, Johns Hopkins University, Baltimore, 2004https://www.ce.jhu.edu/bschafer/cFSM. [8] S. Ádány, Constrained Shell Finite Element Method for Thin-Walled Members, Part 1: Constraints for a Single Band of Finite Elements, Thin-Walled Structures, 2017. [9] S. Ádány, D. Visy, R. Nagy, Constrained Shell Finite Element Method, Part 2: Application to Linear Buckling Analysis of Thin-Walled Members, Thin-Walled Structures, 2017. [10] N. Silvestre, D. Camotim, N.F. Silva, Generalized beam theory revisited: from the kinematical assumptions to the deformation mode determination, Int. J. Struct. Stab. Dyn. 11 (05) (2011) 969–997. [11] R. Gonçalves, M. Ritto-Corrêa, D. Camotim, A new approach to the calculation of cross-section deformation modes in the framework of generalized beam theory, Comput. Mech. 46 (5) (2010) 759–781. [12] N. Silvestre, D. Camotim, Shear deformable generalized beam theory for the analysis of thin-walled composite members, J. Eng. Mech. 139 (8) (2013) 1010–1024. [13] R. Bebiano, R. Gonçalves, D. Camotim, A cross-section analysis procedure to rationalise and automate the performance of GBT-based structural analyses, ThinWalled Struct. 92 (2015) 29–47. [14] R. Gonçalves, P.B. Dinis, D. Camotim, GBT formulation to analyse the first-order and buckling behaviour of thin-walled members with arbitrary cross-sections, ThinWalled Struct. 47 (5) (2009) 583–600. [15] N. Silvestre, D. Camotim, First-order generalised beam theory for arbitrary orthotropic materials, Thin-Walled Struct. 40 (9) (2002) 755–789. [16] N. Silvestre, D. Camotim, Second-order generalised beam theory for arbitrary orthotropic materials, Thin-Walled Struct. 40 (9) (2002) 791–820. [17] D. Camotim, C. Basaglia, N. Silvestre, GBT buckling analysis of thin-walled steel frames: a state-of-the-art report, Thin-Walled Struct. 48 (10–11) (2010) 726–743. [18] C. Basaglia, D. Camotim, N. Silvestre, Global buckling analysis of plane and space thin-walled frames in the context of GBT, Thin-Walled Struct. 46 (1) (2008) 79–101. [19] Z. Li, J.C.B. Abreu, J. Leng, S. Adany, B.W. Schafer, Review: constrained finite strip method developments and applications in cold-formed steel design, Thin-Walled Struct. 81 (2014) 2–18. [20] S. Ádány, B.W. Schafer, A full modal decomposition of thin-walled, single-branched open cross-section members via the constrained finite strip method, J. Constr. Steel Res. 64 (1) (2008) 12–29. [21] N.D.M. Djelil, M. Matallah, M. Djafour, Constrained spline Finite Strip Method for thin-walled members with open and closed cross-sections, Thin-Walled Struct. 132 (2018) 302–315. [22] N. Djafour, M. Djafour, A. Megnounif, M. Matallah, D. Zendagui, A constrained

[28]

[29]

[30]

[31]

[32]

[33]

[34]

[35] [36] [37] [38] [39] [40]

[41]

[42]

[43]

[44] [45]

[46]

14

finite strip method for prismatic members with branches and/or closed parts, ThinWalled Struct. 61 (2012) 42–48. M. Nedelcu, GBT-based buckling mode decomposition from finite element analysis of thin-walled members, Thin-Walled Struct. 54 (2012) 156–163. M. Nedelcu, Buckling mode identification of perforated thin-walled members by using GBT and shell FEA, Thin-Walled Struct. 82 (2014) 67–81. S. Ádány, A.L. Joó, B.W. Schafer, Buckling mode identification of thin-walled members by using cFSM base functions, Thin-Walled Struct. 48 (10–11) (2010) 806–817. Z. Li, S. Ádány, B.W. Schafer, Modal identification for shell finite element models of thin-walled members in nonlinear collapse analysis, Thin-Walled Struct. 67 (2013) 15–24. M. Casafont, F. Marimon, M.M. Pastor, Calculation of pure distortional elastic buckling loads of members subjected to compression via the finite element method, Thin-Walled Struct. 47 (6–7) (2009) 701–729. M. Casafont, F. Marimon, M. Pastor, M. Ferrer, Linear buckling analysis of thinwalled members combining the generalised beam theory and the finite element method, Comput. Struct. 89 (21–22) (2011) 1982–2000. R. Bebiano, D. Camotim, R. Goncalves, GBTUL 2.0-A second-generation code for the GBT-based buckling and vibration analysis of thin-walled members, Thin-Walled Struct. 124 (2018) 235–257. B.W. Schafer, CUFSM 5.01—Elastic Buckling Analysis of Thin-Walled Members by the Finite Strip Method and Constrained Finite Strip Method for General End Boundary Conditions, Department of Civil Engineering, Johns Hopkins University, 2018, http://www.ce.jhu.edu/bschafer/cufsm/. M. Khezri, K.J.R. Rasmussen, An energy-based approach to buckling modal decomposition of thin-walled members with arbitrary cross sections, Part 1: Derivation, Thin-Walled Struct. 138 (2019) 496–517. M. Khezri, K.J.R. Rasmussen, An energy-based approach to buckling modal decomposition of thin-walled members with arbitrary cross-sections, Part 2: modified global torsion modes, examples, Thin-Walled Struct. 138 (2019) 518–531. G. Piccardo, G. Ranzi, A. Luongo, A direct approach for the evaluation of the conventional modes within the GBT formulation, Thin-Walled Struct. 74 (2014) 133–145. S. Jin, D. Gan, H. Chen, R. Cheng, X. Zhou, A force-based method for identifying the deformation modes of thin-walled members, Thin-Walled Struct. 129 (2018) 473–487. J. Becque, X. Li, B. Davison, Modal decomposition of coupled instabilities: the method of the equivalent nodal forces, Thin-Walled Struct. 143 (2019) 106229. J. Becque, A new approach to modal decomposition of buckled shapes, Structure 4 (2015) 2–12. ANSYS, Inc, ANSYS Release 17.0 Documentation, ANSYS, Inc., 2015. N. Silvestre, D. Camotim, On the mechanics of distortion in thin-walled open sections, Thin-Walled Struct. 48 (7) (2010) 469–481. SSMA, SSMA Product Technical Guide, Steel Stud Manufacturers Association2015.. M. Abbasi, M. Khezri, K.J.R. Rasmussen, B.W. Schafer, Elastic buckling analysis of cold-formed steel built-up sections with discrete fasteners using the compound strip method, Thin-Walled Struct. 124 (2018) 58–71. Z. Li, B.W. Schafer, Buckling analysis of cold-formed steel members with general boundary conditions using CUFSM: conventional and constrained finite strip methods, 20th International Specialty Conference on Cold-Formed Steel Structures Recent Research and Developments in Cold-Formed Steel Design and Construction, University of Missouri-Rolla, 2010, pp. 17–31. Z. Li, B. Schafer, Finite Strip Stability Solutions for General Boundary Conditions and the Extension of the Constrained Finite Strip Method, Trends in Civil and Structural Engineering Computing, Saxe-Coburg Publications, Stirlingshire, UK., 2009, pp. 103–130. Z. Li, M.T. Hanna, S. Adany, B.W. Schafer, Impact of basis, orthogonalization, and normalization on the constrained Finite Strip Method for stability solutions of open thin-walled members, Thin-Walled Struct. 49 (9) (2011) 1108–1122. B.W. Schafer, Review: the Direct Strength Method of cold-formed steel member design, J. Constr. Steel Res. 64 (7–8) (2008) 766–778. S. Ádány, B.W. Schafer, Generalized constrained finite strip method for thin-walled members with arbitrary cross-section: secondary modes, orthogonality, examples, Thin-Walled Struct. 84 (2014) 123–133. N. Silvestre, D. Camotim, Nonlinear generalized beam theory for cold-formed steel members, Int. J. Struct. Stab. Dyn. 03 (04) (2003) 461–490.