GBT-BASED ELASTIC-PLASTIC ANALYSIS OF COLD-FORMED STEEL MEMBERS

This paper presents a formulation of Generalised Beam Theory (GBT) intended to perform thorough first-order elastic-plastic analyses of thin-walled members subjected to arbitrary deformations and made of an isotropic non-linear material. The J 2 -flow theory is used to model plasticity in conjunction with the Euler-Backward return-mapping algorithm. After presenting the formulation, its application is illustrated by means of the first order analysis of a simply supported Z-section beam made of an elastic-perfectly plastic material (e.g., carbon steel) and acted by a load uniformly distributed along the flanges. The set of GBT-based results comprises the load-deflection curves (equilibrium paths), displacement profiles, stress distributions (diagrams and 3D contours), and deformed shapes (modal amplitude functions and 3D configurations). These results are compared with the ones obtained from shell finite element analyses (SFEA) using A BAQUS . It is seen that the GBT results display a very good agreement with the SFEA values.


INTRODUCTION
The use of high performance materials, such as high-strength and high-ductility metals, has led to an economy in the weight of built structures.In fact, steel structures are progressively more slender and instability phenomena play an important role in their design.Steel members often display thin walls and slender cross-sections, exhibiting high susceptibility to global (flexural and torsional) and local buckling.An accurate prediction of the buckling and collapse behaviour of steel members has been possible by means of the performance of experimental tests.However, they are expensive and their set-up requires a careful preparation of all devices.A much cheaper alternative is the performance of SFEA, using nonlinear constitutive laws and incremental-iterative techniques.However, the SFEA output data is generally difficult to interpret because (i) it is based on stresses and strains, not stress resultants, and (ii) it provides member deformed configurations with mixed global-local modes.Since the use of complex non-linear 2 SFEA, as well as the interpretation of the results obtained, are time-consuming tasks, a very promising alternative to this approach is the use of Generalised Beam Theory (GBT).Currently, GBT is widely recognized by the scientific community as a powerful, versatile, elegant and efficient thin-walled beam theory.Its elegance arises from the modal nature, i.e., the displacement field is obtained as a linear combination of cross-section deformation modes that vary along the member length.GBT has attracted the interest of several researchers, leading to the development of new formulations and applications.The theory has been extensively upgraded at the Technical University of Lisbon [1] and has been applied to (i) different analysis (first order, buckling, vibration, post-buckling) and (ii) distinct materials (steel, steelconcrete, FRP).In these works, the material was always assumed elastic, with no degradation (plasticity) involved.The first physically non-linear GBT formulation was presented by Gonçalves and Camotim [2], in the context of elastic-plastic bifurcation analysis.Recently, Gonçalves and Camotim [3,4] proposed a finite element based on the J 2 -flow plasticity theory in the context of GBT.In parallel, Abambres et al. [5] developed a GBT first order elastic-plastic formulation, based on the J 2 -flow plasticity theory, to analyse members undergoing only global (axial, flexural, torsional) deformation.In order to overcome this limitation, the work presented in this paper extends this formulation to account for arbitrary (local and global) deformations.Although J 2 -flow plasticity theory was already adopted by Gonçalves and Camotim [3,4], the GBT formulation presented herein has novel features, mainly because it is based on a different set of deformation modes, obtained through the GBT cross-section analysis developed by Silva et al. [8].In order to illustrate the application of the proposed formulation, a simply supported cold-formed steel Zsection beam is analysed and the GBT results are thoroughly discussed and compared with those yielded by ABAQUS [6] SFEA.

BRIEF OVERVIEW OF THE GBT FORMULATION
Due to space limitations, only a brief overview of the GBT formulation is presented in this section.Detailed information can be found in Abambres et al. [5,7].The GBT analysis of a structural member consists of two main steps: (i) a cross-section analysis and (ii) a member analysis.The cross-section analysis comprises the determination of the deformation modes and corresponding displacement profiles (u k (s), v k (s), w k (s)) along the cross-section mid-line (coordinate s -see Fig. 1).The member analysis comprises the determination of the mode amplitude functions  k (x) that characterise the variation of the mode displacement profiles along the member axis x, i.e., the GBT displacement field corresponds to a linear combination of cross-section deformation modes (u=u k (s) k,x (x), v=v k (s) k (x), w=w k (s) k (x)).As mentioned earlier, the cross-section analysis adopted in this work is based on the approach developed by Silva et al. [8], wherein four deformation mode families (transverse extension, warping shear, local and global) are obtained via the solution of a sequence of eigenvalue problems.All the existing GBT approaches consider four degrees-of-freedom (dof) per node, corresponding to three displacements and one transverse flexural rotation.Besides these four dof, the proposed GBT formulation adopts, for the first time, an additional dof: a rotation orthogonal to each wall mid-plane, designated as warping rotation.It is considered in each node  and wall direction, making it possible to approach the warping profile u(s) by means of piecewise Hermite cubic polynomials, instead of the usual piecewise linear functions [8].In what follows, the derivation of the GBT incremental equilibrium equations is briefly described.
It is well known that the determination of a non-linear equilibrium configuration requires the use of an incremental-iterative strategy.Thus, we consider the GBT incremental formulation, defined by int where f int is the internal force vector, λ is the load parameter and K tan is the tangent stiffness matrix at the j th equilibrium configuration.Furthermore, the component i-k of the tangent stiffness matrix is a 4×4 sub-matrix related to modes i and k [7].In addition, it is worth pointing out that the FE approximation is based on Hermite cubic polynomials (see [5]) and that the expressions presented in [7] are not valid for axial and warping shear modes i .Finally, it is also necessary to define the internal force vector at an arbitrary equilibrium configuration.For all the modes displaying in-plane displacements, the component k of the elementary internal force vector at the j th equilibrium configuration (4×1 sub-vector) is obtained from , where  mn is a generic 2 nd Piolla-Kirchhoff stress component, E mn are sectional deformations (not strains) and Ψ H are Hermite cubic polynomials.Concerning the plasticity model, the J 2 -flow theory with the associated flow rule was adopted [9].In order to compute the tangent stiffness matrix, the consistent elastic-plastic constitutive matrix D ep* was used [7], since it is known that it reduces the computational cost of the incremental-iterative procedure, when compared to the conventional constitutive matrix [9].In order to compute the internal force vector (Eqs.( 2)-( 3)), one must be able to update the stresses and the plastic proportional factor at any equilibrium configuration.For that purpose, a robust return-mapping scheme should be used for those points (configurations) undergoing plastic flow during a deformation increment.In this work, the Euler-backward method based on a Newton-Raphson scheme was considered [9].The GBT formulation was implemented using MATLAB [10].

ILLUSTRATIVE EXAMPLE
Now, an illustrative example is presented and the GBT results are compared with those yielded by ABAQUS SFEA [6].We consider the simply supported Z-section beam depicted in Fig. 2 which has length L=1200 mm and the cross-section (mid-line) dimensions shown in Fig. 2(a), and is acted by a uniformly distributed load p = 0.1λ N/mm 2 (Fig. 2(b)) applied at both flanges (Fig. 2(a)).The coldformed steel is modeled as an isotropic linear elastic-perfectly plastic material with Young's modulus i The components i-k due to those modes are obtained by replacing the corresponding vectors Ψ (Hermite shape functions) and their derivatives by vectors with one order of differentiation below.
E=200000N/mm 2 , Poisson's ratio ν=0.3 and uniaxial nominal yield stress σ y =450N/mm 2 .A total of 9 nodes were used to discretise the Z-section (see Fig. 3), involving 7×5+2×6=47 dof -recall that the 2 nodes located at the Z-section corners have 6 dof, corresponding to 3 displacements, 1 flexural rotation and 2 warping rotations (each of them orthogonal to a converging wall).The eigenvalue problem (cross-section analysis) provides 47 deformation modes: 4 global modes, 16 local modes, 19 warping shear modes and 8 transverse extension modes.Due to the minute contribution of 10 deformation modes, only the remaining 37 modes were included in the member analysis  the most relevant are represented in Fig. 3. Regarding the FE models, (i) 15 beam FEs were used in the GBT analysis and (ii) 94 "transverse rows" of 14 shell FEs were considered to obtain the ABAQUS results.
Concerning the numerical integration, 5 Gauss points were used by both the GBT FE model, in all directions (s, z, x) and the SFE model, in the through-thickness direction.that the most relevant modes are mode 2 (major axis bending), with a participation between 12.53-13.28%,mode 3 (minor axis bending), with 73.30-76.46%participation, and mode 5 (local), with 8.24-9.25% participation -note that the participations of these modes change slightly along the whole equilibrium path.It should be noticed that mode 7, which has a negligible participation in the elastic regime, increases its contribution in the elastic-plastic and plastic phases to 3.83%.The ABAQUS and GBT beam deformed shapes at collapse (configuration C), shown in Fig. 4(b), are very similar -note that the deformed shape on the right side of Fig. 4(b) is not a SFE representation, but a 3D view obtained from a beam theory that includes local deformations (GBT).Despite the load being vertical, it is interesting to highlight the huge lateral (horizontal) displacement at failure, mostly due to the participation of mode 3 (minor-axis bending).Notice also the local deformation at mid-span, involving flange rotations in opposite senses and web transverse bending, mostly due to the contribution of modes 5 and 7.   3(d))) were also included in the analysis (they play an important role in the beam behaviour, especially in the plastic regime), their amplitudes are negligible when compared to those of global, local and shear modes.In configuration A (elastic -Fig.5(a)), (i) the warping shear modes 21 and 23 vary linearly along the beam axis, in accordance with classical beam theory (note that mode 21 contributes to the web shear deformations and mode 23 contributes to both the web and flange shear deformation  see Fig. 3(c)), (ii) the global modes 2 and 3 vary cubically along the beam axis, also in accordance with classical beam theory, and (iii) the local modes 5 and 7 exhibit a very smooth variation: they vary considerably near the supports (0≤x≤150mm; 1050≤x≤1200mm) but remain almost uniform in a large fraction of the beam span (150≤x≤1050mm)  note that the variation of mode 3 is larger than that of mode 2, due to the lower second moment of area.In configuration C (plastic -Fig.6 beam axis, but only in the elastic and elastic-plastic zones (0≤x≤400mm; 800≤x≤1200mm), showing a very pronounced non-linear variation with x where plasticity is spreading (400≤x≤800mm), (ii) the global modes 2 and 3 vary almost linearly in the elastic and elastic-plastic zones (0≤x≤400mm; 800≤x≤1200mm), but show a high non-linear variation in the mid-span zone (400≤x≤800mm), where plasticity controls, and (iii) the participations of local modes 5 and 7 decrease in the elastic and elastic-plastic zones (0≤x≤400mm; 800≤x≤1200mm), but increase in the mid-span plastic zone (400≤x≤800mm).This shows that local modes influence the beam plastic resistance, as they are associated with very strong web transverse bending, which erodes the beam strength in the mid-span zone, as can be observed in Fig. 4(b).Since the warping displacements are given by u=u k (s) k,x (x), Fig. 7 shows the variation of the warping shear modes (already depicted in Fig. 6(a)), as well as the derivatives  k,x (x) of global modes 2 and 3, for configuration C (plastic).Note that the warping displacements remain almost unchanged in the elastic and elastic-plastic zones (0≤x≤400mm; 800≤x≤1200mm) but vary a lot (even change sign) in the plastic midspan zone (400≤x≤800mm)  they are null at the mid-span section (x=600mm).For the three equilibrium configurations (A, B, C), Fig. 8(a) depicts the variation of the vertical displacement δ of the top flange edge node (Fig. 2(a)) along the beam axis (0≤x≤1200mm).While the GBT and ABAQUS displacement profiles are nearly coincident at A and C, a maximum difference of 3% occurs at B. Note that the curve (x) concerning configuration C is almost straight for 0≤x≤400mm and 800≤x≤1200mm but is quite non-linear in the mid-span zone (400≤x≤800mm).This shows that the former zones are still in the elastic or elastic-plastic regimes, while the latter is fully yielded (forms a "plastic hinge").Fig. 8(b) displays the Mises stress diagrams at mid-span for the configurations A, B and C. Once again, the GBT and ABAQUS results exhibit a quite satisfactory agreement.In configuration A (elastic), the maximum stress takes place at the web-flange node, while nearly null values occur at the web and flange mid-points.In configuration B (elastic-plastic), the minimum stress takes place at the web central and flange edge zones, while it equals the yield stress in the remaining Z-section zones.In configuration C (plastic), the stress equals the yield stress in practically the whole Z-section (except in part of the flanges, for the ABAQUS diagram)., show that the neutral axis passes through the Z-section centroid and crosses both flanges -note that the flange tips and flange internal zones yield under axial stresses with opposite signs.It is also clear that the neutral axis location remains practically unaltered as the beam evolves from configuration A (elastic) to configuration C (plastic).The shear and Mises stress contours at configuration C (plastic), presented in Figs. 10 and 11, provide further evidence of the similarity between the ABAQUS and GBT results.The sole exception is the small fluctuation of the GBT shear stress magnitude near mid-span, which disappears if a more refined FE is adopted at the mid-span zone.Fig. 11 confirms the previous assertion that the mid-span zone (400≤x≤800mm) is fully yielded, while the outer zones remain elastic (0≤x≤200mm and 1000≤x≤1200mm) or elastic-plastic (200≤x≤400mm and 800≤x≤1000mm).

CONCLUSIONS
A GBT elastic-plastic formulation was employed to perform a first order elastic-plastic analysis of a cold-formed simply supported Z-section steel beam subjected to a uniformly distributed load.The results presented consisted of (i) load-deflection curves (equilibrium paths), (ii) displacement longitudinal profiles, (iii) axial, shear and Mises stress diagrams, including 3D contours, and (iv) 3D deformed configurations.In order to highlight the modal nature of GBT, the variation of modal participation along the equilibrium path was shown, as well as the shapes of the modal amplitude functions  for this purpose, GBT results concerning three distinct equilibrium configurations (elastic, elastic-plastic and plastic) were presented.All these results, found to be in very good agreement with those provided by ABAQUS (the best correlation was obtained for the equilibrium paths and displacement profiles  differences always below 3%), were discussed in great detail and revealed the potential of the proposed approach to provide in-depth insight on the spread of plasticity and formation of plastic mechanisms.Finally, it is important to highlight that, in the illustrative example shown, the number of d.o.f required by the ABAQUS SFEA was 7.4 times higher than that involved in using the GBT model, which attests the computational efficiency of GBT analyses.

Fig. 4 (
Fig.4(a) presents the equilibrium paths λ(δ) obtained from the GBT (circles) and ABAQUS (triangles) analyses, where δ is the vertical displacement at the mid-span top flange edge node (Fig.2(a)) (x=600 mm -Fig.2(b)).In order to compare the GBT and ABAQUS results, three equilibrium configurations were selected: A (elastic), B (elastic-plastic) and C (plastic).It is seen that the GBT and SFE (ABAQUS) curves compare very well, with a maximum difference of 1.5%.In particular, the load parameter values at the 3 equilibrium configurations A, B and C (see Fig.4(a)) are (i) λ GBT =0.39 and λ SFE =0.39 (A), (ii) λ GBT =3.38 and λ SFE =3.36 (B), and (iii) λ GBT =3.63 and λ SFE =3.59 (C).Using GBT, the yield load is p y 0.233 N/mm 2 (2.33) and the ultimate load is p u =0.363 N/mm 2 (=3.63),corresponding to a shape factor S=p u /p y =1.56.The GBT modal participation diagram is displayed in the bottom of Fig.4(a) and shows

Figs. 5
Figs. 5(a)-(b) and 6(a)-(b) show the variation of the modal amplitude functions along the beam length: while Figs.5(a) and 6(a) provide the amplitude variation (ζ k,x (x)) of the warping shear modes 21 and 23 (see Fig. 3(c)), Figs.5(b) and 6(b) depict the variation (ζ k (x)) of the global (2 and 3 -see Fig. 3(a)) and local (5 and 7 -see Fig. 3(b)) modes.Although transverse extension modes (see Fig.3(d))) were also included in the analysis (they play an important role in the beam behaviour, especially in the plastic regime), their amplitudes are negligible when compared to those of global, local and shear modes.In configuration A (elastic -Fig.5(a)),(i) the warping shear modes 21 and 23 vary linearly along the beam axis, in accordance with classical beam theory (note that mode 21 contributes to the web shear deformations and mode 23 contributes to both the web and flange shear deformation  see Fig.3(c)), (ii) the global modes 2 and 3 vary cubically along the beam axis, also in accordance with classical beam theory, and (iii) the local modes 5 and 7 exhibit a very smooth variation: they vary considerably near the supports (0≤x≤150mm; 1050≤x≤1200mm) but remain almost uniform in a large fraction of the beam span (150≤x≤1050mm)  note that the variation of mode 3 is larger than that of mode 2, due to the lower second moment of area.In configuration C (plastic -Fig.6(a)),(i) the warping shear modes 21 and 23 still vary linearly along the Figs. 5(a)-(b) and 6(a)-(b) show the variation of the modal amplitude functions along the beam length: while Figs.5(a) and 6(a) provide the amplitude variation (ζ k,x (x)) of the warping shear modes 21 and 23 (see Fig. 3(c)), Figs.5(b) and 6(b) depict the variation (ζ k (x)) of the global (2 and 3 -see Fig. 3(a)) and local (5 and 7 -see Fig. 3(b)) modes.Although transverse extension modes (see Fig.3(d))) were also included in the analysis (they play an important role in the beam behaviour, especially in the plastic regime), their amplitudes are negligible when compared to those of global, local and shear modes.In configuration A (elastic -Fig.5(a)),(i) the warping shear modes 21 and 23 vary linearly along the beam axis, in accordance with classical beam theory (note that mode 21 contributes to the web shear deformations and mode 23 contributes to both the web and flange shear deformation  see Fig.3(c)), (ii) the global modes 2 and 3 vary cubically along the beam axis, also in accordance with classical beam theory, and (iii) the local modes 5 and 7 exhibit a very smooth variation: they vary considerably near the supports (0≤x≤150mm; 1050≤x≤1200mm) but remain almost uniform in a large fraction of the beam span (150≤x≤1050mm)  note that the variation of mode 3 is larger than that of mode 2, due to the lower second moment of area.In configuration C (plastic -Fig.6(a)),(i) the warping shear modes 21 and 23 still vary linearly along the

Fig. 5 .
Fig. 5. Variation of modal amplitude functions at configuration A: (a) warping shear and (b) global and local modes.

Fig. 6 .
Fig. 6.Variation of modal amplitude functions at configuration C: (a) warping shear and (b) global and local modes.

Fig. 7 .
Fig. 7. Derivative of the amplitude functions ζ k (x) of the global and warping shear modes (configuration C).

Fig. 8 .
Fig. 8. (a) Longitudinal variation of the vertical displacement of the flange edge node and (b) Mises stress (σ Mises ) diagram for configurations A, B and C (mid-span section x=596 mm).

Figs. 9
Figs. 9(a)-(b) provide the axial stress diagrams at the mid-span section for configurations A and B (Fig.9(a)), and C (Fig.9(b)), and show that there is a quite satisfactory resemblance between the ABAQUS and GBT results.From the observation of these figures, it is concluded that plasticity spreads from the web-