International Journal of Solids and Structures 48 (2011) 1894–1905

Contents lists available at ScienceDirect

International Journal of Solids and Structures journal homepage: www.elsevier.com/locate/ijsolstr

Static behavior and bifurcation of a monosymmetric open cross-section thin-walled beam: Numerical and experimental analysis Angelo Di Egidio a,⇑, Fabrizio Vestroni b a b

Department of Structural, Hydraulic and Geotechnical Engineering, University of L’Aquila, L’Aquila, Italy Department of Structural and Geotechnical Engineering, University of Rome ‘‘La Sapienza’’, Rome, Italy

a r t i c l e

i n f o

Article history: Received 16 December 2010 Received in revised form 14 February 2011 Available online 3 March 2011 Keywords: Open cross section beam Nonlinear warping Torsional elongation Experimental investigation Finite element model Flexural–torsional instability

a b s t r a c t The aim of the paper is the numerical and experimental validation of a previously developed nonlinear one-dimensional model of inextensional, shear undeformable, thin-walled beam with an open cross-section. Nonlinear in-plane and out-of-plane warping and torsional elongation effects are included in the model. To better understand the role of these new contributions a beam with a section with one symmetry axis, undergoing moderately large flexural curvatures and large torsional curvature is taken into account. To obtain a section of a cantilever beam for which the torsional curvature is expected to prevail with respect to the flexural ones, a preliminary study is performed. The attention is focused on the response to static forces and on the stability of the equilibrium branches. Analytical results are compared with results of two different nonlinear finite element models and mainly with experimental results to confirm the validity of the analytical model. Interesting results are obtained for the critical values of the flexural–torsional instability loads. Ó 2011 Elsevier Ltd. All rights reserved.

1. Introduction The description of the nonlinear behavior of beams in 3D has been the subject of a certain number of papers in the last decades. A one-dimensional polar model of a compact beam was initially studied in Crespo da Silva and Glynn (1978); the torsional component was statically condensed and warping was neglected. A threedimensional beam model was developed in Rosen and Friedmann (1979) for a compact cross-section beam by also taking into account the linear warping. In Crespo da Silva (1988, 1991) and Crespo da Silva and Zaretzky (1994) flexural–torsional–extensional couplings in the motion of a cantilever beam was considered, limiting the model to linear warping. Even if the linear warping contribution is considered, attention is generally paid to the case of compact cross-sections. In Luongo et al. (1989) a non-compact cross-section beam, having close bending and torsional frequencies was studied but due to difficulties in tackling nonlinear warping, it was preferred to completely neglect it. An approach based on the extension of the Vlasov theory (Vlasov, 1961) to the nonlinear field is found in Ghobarah and Tso (1971), Mollmann (1982) and Ascione and Feo (1995); however, due to the complexity of the problem, several simplifying assumptions have been used, i.e. flexural curvatures up-to linear terms are considered. In recent years sim-

⇑ Corresponding author. E-mail addresses: [email protected] (A. Di Egidio), [email protected] (F. Vestroni). 0020-7683/$ - see front matter Ó 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.ijsolstr.2011.02.029

ilar problems have been studied in a numerical way by using the model developed in Schulz and Filippou (1998) where a non-uniform warping theory of bars employing two-dimensional St. Venant warping functions and one-dimensional independent warping parameters has been developed. In Sapountzakis and Mokos (2003) and Sapountzakis and Tsiatas (2007, 2008), by using boundary value problems through FEM, warping and shear stresses of bars under arbitrary loading have been evaluated, while static and dynamic analyses accounting of nonlinear torsional have been performed. More recently in Di Egidio et al. (2003, Part I and II) the description of the mechanical behavior of beams with open cross-sections is dealt with by using a nonlinear beam model developed in a rigorous manner starting from an internally constrained threedimensional continuum, in which torsional and flexural curvatures of the same order of magnitude are considered. The warping is obtained by extending the Vlasov theory (Vlasov, 1961) to the nonlinear field. The effects of the torsional curvature on the elongation of the longitudinal fibers and the nonlinear warping of the section are considered. The model obtained is very complex and the equilibrium equations are derived by means of symbolic manipulation tools missing any possibility of recognizing the mechanical meaning of the different numerous terms. Within a class of open cross section beams with certain geometrical and mechanical characteristics, it can be shown that torsional curvature is greater than the flexural ones. A notable simplification is then obtained in kinematical relations with respect to the model developed in Di Egidio et al. (2003, Part I and II), as partially presented in Vestroni et al. (2006).

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

In this paper three simplified static equations are derived describing the behavior of inextensional and shear undeformable nonlinear 3D beam with a section that presents one symmetry axis. The new terms due to nonlinear warping and torsional elongation are now well recognizable with respect to the equations obtained in Di Egidio et al. (2003, Part I and II). A preliminary study is developed to determine the geometrical and mechanical characteristics of a cantilever beam with a channel section, subject to a force acting on the centroid of the free-end section, for which the beam undergoes deformations with the torsional curvature greater than the flexural ones. A Galerkin discretization is performed by using different shape functions. Special care is devoted to the choice of these discretizing functions. As for the study of moving masses on a beam (Biondi and Muscolino, 2005) and for the study of the linear stability of a beam with a concentrated dashpot (Luongo and D’Annibale, submitted for publication), static functions are used. These functions make it possible to better represent the effects of a concentrated force, since they make it easier to describe the discontinuity of the shear force. The attention is focused on the response to static forces and on the stability of the equilibrium branches. Analytical results are compared with results of four different nonlinear FE models and mainly with experimental results. Numerical and experimental investigations are carried out to confirm the importance of the new nonlinear contributions due to warping and torsional elongation and to validate the model developed in Di Egidio et al. (2003, Part I). The good performance of the simplified model here proposed, also permits the selection of the prevailing terms among the many new nonlinear ones of model proposed in Di Egidio et al. (2003, Part I) and the better understanding of their kinematical nature.

2. Analytical model An initially straight thin-walled beam with an open cross-section, arbitrarily restrained at the ends, is considered. A model is developed which is accurate enough to describe the behavior of a slender monosymmetric channel section beam undergoing large displacements, but at the same time, rather simple to make possible a mechanical interpretation of the main contributions in the equations. 2.1. Kinematics The following hypotheses are assumed: (a) the beam cross-section is rigid and remains orthogonal to the centroid axis in the deformed configuration (shear indeformable beam); (b) a non-rigid displacement field is superimposed on the previous one, having components both normal and tangential to the cross-section in the deformed configuration (warping); (c) the shear strains on the middle surface of the thin-walled beam identically vanish (Vlasov condition); moreover, the extensional and shear strains of the cross-section plane also vanish (indeformability of the section in its own plane); (d) the beam is axially inextensible; (e) the torsional curvature l3 is prevailing with respect to the flexural curvatures l1, l2. A reference frame (O x1 x2 x3) is introduced, where x1 and x2 are section principal axes, x3 contains the centroid axis and O is placed on the centroid G of an end cross-section. A unit base vector b = {b1, b2, b3}T, solid with the (not warped) section in the deformed configuration is considered, with b3 tangent to the centroid-axis. Let denote by b = {b1, b2, b3}T the triad solid with the section in the undeformed configuration, oriented like the xi-axes (Fig.1). The displacement vector of the generic point P = (x1, x2,

1895

Fig. 1. Beam section before and after deformation, and unit vector triads.

x3) can be expressed as the sum of a rigid and a non-rigid displacement namely:

uP ¼ uC þ ðR  IÞðx  xC Þ þ R/

ð1Þ

where R is the rotation matrix, I the identity matrix, u = {u1 u2 u3}T the displacements of the shear center C  ðx1C ; x2C Þ; x ¼ fx1 x2 0gT and / = {/1 /2 /3}T the warping functions. To make the kinematical description unique, the warping vector must describe neither a translation nor a rotation. This requirement is satisfied if the following orthogonality conditions hold:

Z /2 dA ¼ 0; /3 dA ¼ 0 A A Z ZA /3 x1 dA ¼ 0; /3 x2 dA ¼ 0; A ZA ½/2 ðx1  x1C Þ  /1 ðx2  x2C ÞdA ¼ 0

Z

/1 dA ¼ 0;

Z

ð2Þ

A

Eq. (2)1–3 ensure / is not a translation; Eq. (2)4,5 assure / is not a rotation for any /3 – 0; Eq. (2)6 prevents / from being a purely torsional rotation. The orthogonal rotation matrix R describes the position of the base b with respect to the base b (b = Rb). The displacements field (1) is described by six functions of the abscissa z :¼ x3 and of the time t, i.e. ui(z, t) and 0i(z, t) and by the three warping functions /i(x1, x2, z, t) (i = 1, 2, 3). The antisymmetric curvature matrix C, referred to the undeformed base b is C = RTR0 where ()0 = @/oz. The angular velocity matrix W, referred to the undeformed base b, is also an antisymmetric matrix and it is given by a similar relation W ¼ RT R_ where the dot denotes timedifferentiation. The Green–Lagrange strain matrix E = [eij] is assumed as the deformation measure. Under the hypothesis of shear indeformability it is possible to express the strains eij as functions of the curvatures li and of the elongation eC of the shear center axis, i.e. of the generalized strain measures of the one-dimensional polar beam model, in addition to the derivatives of the warping functions. The attention is focused on cases for which the beam undergoes large torsional curvature and moderately large flexural curvatures; this assumption will make it possible to analyze the flexural–torsional instability of the beam and to compare results with those known in the technical literature. The ordering l3 = O(e), li = O(e3/2) (i = 1, 2) is then introduced with e a small parameter. Assuming eC and /i (i = 1, 2, 3) of order e, the Green–Lagrangian strains, corrected up-to the e2-order, read:

1896

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

e11 ¼ /1;1 þ 1=2/23;1 e22 ¼ /2;2 þ 1=2/23;2 c12 ¼ /1;2 þ /2;1 þ /3;1 /3;2 c13 ¼ l3 ðx2  x2C Þ þ /3;1 þ /1;3 c23 ¼ l3 ðx1  x1C Þ þ /3;2 þ /2;3 e33 ¼ eC þ l1 ðx2  x2C Þ  l2 ðx1  x1C Þ þ /3;3 þ 1=2l23 ½ðx1  x1C Þ2 þ ðx2  x2C Þ2 

Z

ð3Þ

02 where eC ¼ u03 þ 1=2ðu02 1 þ u2 Þ (see Di Egidio et al. (2003, Part I) for more details) and cij = 2eij. The term proportional to the squared torsional curvature l3 appearing in the longitudinal strain e33 represents the nonlinear elongation of the longitudinal fibers of the beam due to the torsional rotation of the section. Finally, the inplane strains only depend on warping, since the cross-section is rigid in bending and torsion.

2.2. Nonlinear warping In order to obtain a one-dimensional model, the dependence of the warping functions /i(x1, x2, z, t) on the transversal coordinates x1 and x2, should be determined. With the aim to evaluate the inplane warping components /1 and /2, the in-plane strains e11, e22 and c12 (Eq. (3)1–3) are required to vanish on the whole section. In general no function having this property exists, since the plane strain problem is overdetermined; however, /1 and /2 exist if /3 is taken constant along the thickness, as usually done for thin-walled beams, and can themselves be considered constant along the thickness. This fact makes the warping components dependent only on the local abscissa c (see Fig. 1) running along the middle line of the thin-walled section. Finally to evaluate the warping functions /k (k = 1, 2, 3) the strains ecc, enn, ccn and c3c are first evaluated, where n is the inward normal to the middle line at P. By enforcing the kinematical conditions ecc = 0, enn = 0, czc = 0 and cc3 = 0, and by imposing the conditions given by Eq. (2), the following final form of the warping is obtained

/1 ¼ b1 l23 ;

/2 ¼ b2 l23 ;

/3 ¼ a1 l3 þ b3 l03 l3

ð4Þ

2

where e -order terms are only retained, a1 is the Vlasov sectorial area and the new functions bi are defined in Appendix A. By replacing Eq. (4) in Eq. (3)6, the only non-vanishing strain e33 is determined:

e33 ¼ eG þ l1 x2  l2 x1 þ l

0 3

1 a1 þ l23 s2 þ ðl3 l03 Þ0 b3 2

ð5Þ

where eG ¼ eC  l1 x2C þ l2 x1C is the longitudinal strain of the centroid axis and s2 is the square of the distance between the shear center of the section and the generic point P on the middle line. The first four terms of Eq. (5) are formally equal to those of the Vlasov linear theory, while the flexural and torsional curvatures are nonlinear. The fifth term describes the elongation due to torsion, and the remaining term accounts for nonlinear warping. 2.3. Internal constraints The shear indeformability condition makes it possible to express the flexural rotations #1 and #2 in terms of the spatial derivatives of the displacements ui, thus reducing the number of independent displacement variables. The following well-known relationships are drawn (Crespo da Silva and Glynn, 1978):

u02 tg #1 ¼  ; 1 þ u03

u01 tg #2 ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 0 2 u02 1 þ ð1 þ u3 Þ

e33 dA ¼ 0

ð7Þ

A

ð6Þ

The beam is also assumed to be inextensible. The inextensibility condition is introduced by requiring the mean value on the crosssection of the longitudinal strain (5) to vanish, i.e.:

From Eq. (7) it is possible to obtain u03 ¼ f ðu1 ; u2 ; #3 ; tÞ when shear indeformability conditions (6) are used to express #1 and 02 as functions of the other displacement components. Finally, Eq. (5) can be written as e33 = ef + et + e/, where the subscript f stands for flexural, t stands for torsional elongation and / for warping and

1 2

ef ¼ l1 x2  l2 x1 ; et ¼ l23 ðs2  q2C Þ; e/ ¼ l03 a1 þ ðl3 l03 Þ0 b3 ð8Þ in which qC is the cross-section polar radius of inertia with respect to the shear center C. Expanding in series the nonlinear curvatures, taking into account the internal constraints (6) and (7) and remembering the ordering of the curvatures, for sections having x1 as symmetry axis ðx2C ¼ 0Þ, the flexural and torsional curvatures lead to the following relationships:

l1 ¼ u002 þ u001 #3  x1C ðu001 u2 þ u0001 u02 Þ 0 000 l2 ¼ u001 þ u002 #3 þ x1C ðu002 1 þ u1 u1 Þ l3 ¼ #03  u01 u002  x1C ð2u01 u001 u002 þ u001 u0001 u02 Þ

ð9Þ

In the flexural curvatures l1 and l2 the only quadratic terms are retained since it has been considered that they are of order e3/2 while the torsional curvature is of order e. 2.4. Static equations In this section, static equations of the thin-walled beam with the axis x1 as symmetry axis, are obtained. By restricting the analysis to isotropic beams and neglecting the effect of the Poisson ratio, the elastic potential energy per unit length reads:

Z 1 ½Gðc231 þ c232 Þ þ Ee233 dA 2 A Z 1 1 ¼ GJ l23 þ E ½e2f þ e2t þ e2/ þ 2ðef et þ ef e/ þ et e/ ÞdA 2 2 A

V¼

ð10Þ

After the integration on the area of the section, Eq. (10) reads:

V¼

1 GJl23 2 " # 0 02 2ECft l2 l23 1 EI1 l21 þ EI2 l22 ECt l43 EC1 l02 3 þ EC2 ðl3 l3 Þ þ þ þ þ 2 ef et e2f e2t e2/ ð11Þ

where G and E are elastic moduli and GJ is the St. Venant torsional stiffness obtained by assuming that, in Eq. (11), c31 and c32, (which are not zero out of the section middle line) contribute to the St. Venant torsional elastic term only. Many terms vanish after the integration on the area of the section due to symmetry of the section with respect to the x1 -axis, (i.e. efe/, ete/). Below the underline, the deformations from which the terms arise are pointed out to keeping track of the origin of the terms in the next equations. Coefficients EIi (i = 1, 2) are the flexural stiffness and EC1 is the Vlasov warping torsional stiffness. The new coefficients are defined as follows:

C2 ¼ t

Z

b3 ðcÞdc; Ct ¼ 1=2t c Z Cft ¼ 1=2t sðcÞ2 x1 ðcÞdc

Z c

ðsðcÞ2  q2C Þdc; ð12Þ

c

where t is the small thickness of the section. Finally, the three static equations up-to cubic terms finally read:

1897

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905 000 00 002 000 0 0002 0 EI2 f½x1C ðu001 u000 1 þ u1 u2 #3 þ x1C u1 u1 þ x1C u1 u1 Þ 00 00 0 000 003 þ ½u001 þ #3 u002 þ x1C ð3u002 1 þ 2u1 u2 #3 þ u1 u1 þ 2x1C u1 00 000 0 00 0 00 0 002 þ 2x1C u01 u001 u000 1 Þ  ½x1C ðu1 u1 þ u1 u2 #3 þ x1C u1 u1 þ x1C Þ g 00 00 000 0 þ EI1 f½u002 #3 þ u001 #23 þ x1C ðu002 2  2u1 u2 #3  u1 u2 #3 000 0 00 00 0 00 00 0 þ x1C u001 u002 2 þ x1C u1 u2 u2 Þ  ½x1C ðu2 u2  u1 u2 #3 02 000 0 002 00 0 þ x1C u001 u02 u002 þ x1C u000 1 u2 Þ g þ GJf½u1 u2  u2 #3 0 0 0 0 00 0 00 0 0 0 000  x1C ð2u001 u002 #03 þ u000 1 u2 #3 Þ  ½2x1C u1 u2 #3  þ ½x1C u1 u2 #3  g 0 0002 000 00 þ EC1 f½u001 u002 u000 2 þ u1 u2  u2 #3 0000

00 00 00 000 00 0 00 0  x1C ð3u000 1 u2 #3 þ 2u1 u2 #3 þ u1 u2 #3 Þ 00 00 0 00 000 00 00 00 000 0 00 0 000 00 00 þ ½u001 u002 2  u2 #3 þ x1C ðu1 u2 u2  4u1 u2 #3  u1 u2 #3  2u1 u2 #3 Þ 0000

þ ½x1C ðu001 u02 #003 þ 3u01 u002 #003 Þ000  ½x1C u01 u02 #003  g 00 02 00 0 02 000 þ ECft f½#02 3 þ 2x1C u1 #3   ½x1C u1 #3  g ¼ Q 1

ð13Þ

(EC2 terms). The hypothesis of flexural curvatures smaller than the torsional one makes possible to isolate the prevailing new nonlinear effects and permits to keep track of the kinematical effects from which they arise with respect to the results obtained in Di Egidio et al. (2003, Part I). Finally, it is interesting to observe that cubic terms in the previous equations, depending on the coefficient ECft (related to torsional elongation), are generated by quadratic terms of the flexural curvatures. If only linear terms in these curvatures are retained, in Eq. (14) no new terms appear and this model would not be able to correctly evaluate some aspects of the flexural behavior of the beam in the x2 direction. It is possible to assert that the main new nonlinear effects are related to the torsional elongation since it is present in each static equation (Eqs. (13)–(15)). The nonlinear warping contributes only with small quantitative effects. 3. Geometrical characteristics of the beam

00 00 000 00 000 00 0002 0 0 EI1 f½x1C ðu000 1 u2  u1 u1 #3 þ x1C u1 u1 u2 þ x1C u1 u2 Þ 000 0 002 00 00 000 0 00 þ ½u002  #3 u001 þ x1C ð2u001 u002  u002 1 #3 þ u1 u2 þ x1C u1 u2 þ x1C u1 u1 u2 Þ g 00 0 000 0 000 0 0 þ EI2 f½u001 #3 þ u002 #23 þ x1C ðu002 1 #3  u1 u1 #3 Þ g þ GJf½x1C u1 u1 #3  0000

00 0 00 0 0 0 00 0 00 00 000 00 0 þ ½u02 1 u2  u1 #3  2x1C u1 u1 #3  g þ EC1 f½x1C ðu1 u1 #3 þ u1 u1 #3 Þ 00 00 00 0 00 000 002 00 0 000 00 00 þ ½u002 1 u2  u1 #3 þ u1 u1 u2  x1C ð2u1 #3 þ 3u1 u1 #3 Þ 02 00 000 0 00 00 000  ½u01 u001 u002  u01 #003 þ u02 1 u2  2x1C u1 u1 #3  g þ ECft f½#3 #3  g ¼ Q 2

ð14Þ 00 00 002 00 00 000 0 EI1 fu002 1 #3  u1 u2  x1C ðu1 u2 þ u1 u1 u2 Þg

A preliminary study is performed to determine the geometrical characteristics of the beam for which the hypothesis of torsional curvature larger than the flexural ones is verified when a force F acting on the centroid G of the free-end section of a cantilever beam is considered (Fig. 2). It is well known that the linear static problem of an open cross-section beam is governed by three uncoupled differential equations in the u1, u2 and 03 displacement components. However, a force F applied to the centroid G of the section generates a torsional moment Mt = F(xC + xG).

002 002 0 000 00 þ EI2 fu001 u002 þ u002 2 #3 þ x1C ðu1 u2 þ u1 u1 u2 Þg 0 0 þ GJf½#03  u01 u002  x1C ðu01 u001 u003 þ u01 u000 1 u2 Þ g 00 1 f½#3

00 00 000 0 E    x1C ðu002 1 u2 þ u1 u1 u2 0000 00 00 0 00 000 0 0 3u01 u000 1 u2 þ 2u1 u1 u2 þ u1 u1 u2 Þ g 00 002 000 0 0002 0 002 0 002 E 2 f½#3 #3 þ #3 #3  þ 2½#3 þ #03 #003 #000 3   ½#3 #3

þ C þ

u001 u002

þ C

Table 1 Cross-section characteristics.

u01 u000 2

000 000 þ #02 3 #3  g

0 00 0 00 002 0 0 000 0 0 þ ECft fu002 #03 3  2½u1 #3 þ u2 #3 #3 þ x1C ðu1 #3 þ u1 u1 #3 Þ g 00  ECt f6#02 3 #3 g ¼ Q 3

ð15Þ

h = 0.05 m t = 0.002 m xG = 0.00625 m x1C ¼ 0:0156 m I1 = 8.33  108 m4 C = 5.29  1012 m6 C2 = 5.714  1021 m10

Cft =  3.825  1010 m5

with the relevant boundary conditions, omitted for brevity. Quantities Q1, Q2 and Q3 refer to distributed loads that do work with respect to u1, u2 and #3 displacements. Underlined terms are the new terms due to nonlinear warping and torsional elongation. In particular, in Eqs. (13) and (14) the new terms depend only on the coefficient ECft and arise from the interaction between the flexural and the torsional elongation deformations efet (see Eq. (11)). In Eq. (15) also terms depending on nonlinear warping are present

b = 0.025 m l = 1.0 m xC = 0.00938 m A = 0.0002 m2 J = 2.67  1010 m4 I2 = 1.302  108 m4 q2C ¼ 0:000726 m Ctt = 1.265  1011 m6

Table 2 Material characteristics. E = 2.07  1011 N/m2 m = 1.56 kg/m

Fig. 2. (a) Mechanical system under study; (b) ratio between l3 and l1.

G = 8.61  1010 N/m2 IC ¼ 0:00114 kg m

1898

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

Fig. 3. (a) Linear warping function (Vlasov); (b) nonlinear warping function; (c) elongation shape function.

The following nondimensional parameters at = t/h, ab = b/ h, al = h/l are introduced. By fixing h = 0.05 m and choosing at = 0.04, the ratio between the linear torsional curvature l3 ¼ #03 and the linear flexural curvature l1 ¼ u002 is then analyzed, under the force F. In Fig. 2b a surface representing this ratio for the section in z = l/2 is plotted, as a function of the other two nondimensional parameters ab and al. It is possible to note the existence of a region of the parameters in which the ratio is maximized. The surface in Fig. 2b is topologically similar for each value of the parameter at. Taking in mind the experimental investigation, the following values of the parameters are chosen: ab = 0.5, al = 0.05, corresponding to the point labeled with A on the surface. The ratio between the torsional and the flexural curvatures at point A is equal to 17.72. For these particular values of the nondimensional parameters, and referring to a steel beam, the cross-section and the material characteristics are reported, respectively, in Tables 1 and 2 (Appendix B). It is worth noticing that many commercial channel sections, for certain length l of the beam, are positioned along the crest of the surface in Fig. 2b. To better understand the effects of the warping and torsional elongation effects, the shape functions a1(c), b3(c) and ðsðcÞ2  q2C Þ, appearing in Eq. (8), are represented in Fig. 3. It is possible to note that the nonlinear warping function b3(c) is very smaller than the linear one a1(c) (i.e. Vlasov sectorial area). On the contrary, the torsional elongation effects, given by the shape function ðsðcÞ2  q2C Þ; is of the same order of magnitude of Vlasov linear function. However, when the torsional curvature is not negligible,

the deformation related to the nonlinear warping (Eq. (8))3 can become important (see EC2 terms in Eq. (15)). 4. Numerical and experimental model In this section a description of the numerical and the experimental models is presented. Comparison among the first three frequencies of vibration of the beam, obtained by the different models, is performed. 4.1. Finite element models Four FE models have been built by using two different FE environments: ADINAÒ (Automatic Dynamic Incremental Nonlinear Analysis, http://www.adina.com) and AnsysÒ (http://www.ansys.com), and in each FE environment two different models have been taken into account. In the first model (FEM1) the beam has been discretized by using 25 beam elements with open cross-section, accounting only for the linear Vlasov warping. The second model (FEM2) is a more refined model based on bi-dimensional shell elements. The beam has been discretized by using 100 shell elements; the length of the beam has been divided into 25 parts. These finite element models have been used to evaluate the static deflection of the beam accounting for the geometrical nonlinearities and to evaluate the first critical loads in the flexural– torsional instability.

Fig. 4. (a) Instrumented sections and load; (b) points of the instrumented sections where the displacements are measured; (c) particular of the loaded section.

1899

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

Fig. 5. Experimental model: (a) beam (red box) and instrumented sections (blue boxes); (b, c) loaded free-end section (red boxes). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

4.2. Experimental setup An experimental investigation to evaluate the static behavior of the beam is performed. It is conducted with the objective to point out the role of the new nonlinear contributions due to the warping and the torsional elongation. The load is a torsional moment applied to the free end-section which causes great torsional curvature and small flexural ones; in this case the new effects, mainly dependent on the torsional curvature, can be easily outlined (see Eq. (8)2,3). In Fig. 4 the experimental setup is schematically shown. In Fig. 4a the two instrumented sections (1, 2) and the system adopted to load the free-end section by a torsional moment are shown. In each one of the instrumented sections, in the middle point of their sides an instrument able to measure the orthogonal displacement has been used (Fig. 4b). Finally, in Fig. 4c the loaded sections in the undeformed and in the deformed condition are shown to point out the nature of the torsional moment applied at the end of the beam. Since it is generated by two parallel forces and the distance between these change with the configuration of the section, particular attention will be devoted to model this excitation. In other words, this load represents a non-follower torsional moment. In Fig. 5, pictures of the experimental setup are shown. Instrumented sections (Fig. 5a – blue boxes) and loaded free-end section (Fig. 5b and c – red boxes) are better highlighted. 4.3. Comparison among the frequencies of the different models

€ 1 þ EI2 uIV mu 1 ¼ 0 € 2  mx1C #€3 þ EI1 uIV mu 2 ¼ 0 00 € € 2 þ EC1 #IV IC #3  mx1C u 3  GJ#3 ¼ 0 with the relevant homogeneous boundary conditions:

u1 ð0Þ ¼ 0;

u01 ð0Þ ¼ 0;

EI2 u001 ðlÞ ¼ 0;

EI2 u000 1 ðlÞ ¼ 0

u2 ð0Þ ¼ 0;

u02 ð0Þ #03 ð0Þ

EI1 u002 ðlÞ ¼ 0; E 1 #03 ðlÞ ¼ 0;

EI1 u000 2 ðlÞ ¼ 0

#3 ð0Þ ¼ 0;

¼ 0; ¼ 0;

0 EC1 #000 3 ðlÞ  GJ#3 ðlÞ ¼ 0

C

ð17Þ Since the Eq. (16)1 is uncoupled from the remaining twos, it represents a standard problem of a planar cantilever beam. With reference to the other two coupled Eq. (16)2,3, a separation of variables is performed:



u2 ðz; tÞ #3 ðz:tÞ



 ¼

   /2 ðzÞ ixt A kz ixt e ¼ e e B /3 ðzÞ

ð18Þ

By introducing Eq. (18) into Eq. (16)2,3, the following eigenvalue problem is obtained:

"

EI1 k4  mx2

mx1C x2

#  A

mx1C x2

ECk4  GJk2  IC x2

B

¼

  0 0

ð19Þ

The vanishing of the determinant in Eq. (19) permits to obtain eight values of the spatial frequency k and eight eigenvectors {A B}T expressed as functions of the unknown temporal frequency x. The following solution is then obtained:

 The eigenproblem of the analytical model must be first solved. The linear homogeneous equations of motion read:

ð16Þ

/2 ðzÞ /3 ðzÞ

 ¼

  8 X Ai ðxÞ ki ðxÞz ki e B i ð xÞ i¼1

ð20Þ

1900

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

Table 3 First three frequencies of the system (AM: analytical model; FEM1: finite element model using one-dimensional beam elements; FEM2: FEM using two-dimensional shell elements; EM: experimental model). AM (Hz)

FEM1 (Hz)

FEM2 (Hz)

EM (Hz)

23.13 40.30 81.29

22.35 41.63 80.47

23.39 41.46 82.48

22.89 40.98 82.01

where ki are unknown coefficients of the linear combination. By taking into account Eqs. (18) and (20) and imposing the boundary conditions (17)4–9 the new eigenvalue problem C(x)k = 0 is obtained, where k = {k1, . . ., k8}T and C(x) is a 8  8 matrix. Values of x that vanish the determinant of C, numerically solved, are the eigenvalues of the coupled problem; eigenvectors k completely define the eigenfunctions (20) of the problem. Frequencies obtained by the analytical model (AM) are compared with those furnished by the two different finite element models (FEM1 – AnsysÒ, FEM2 – AnsysÒ) and with those obtained by the experimental measurement (EM). The modal identification has been conducted measuring in some points the response to the hammer impulse force. The identification procedure has been based on the Frequency Domain Decomposition (Brincker et al., 2001). In Table 3 are reported the first three frequencies of the mechanical system obtained by the different models. Results of the analytical model and of the FE models are in good agreement with the frequencies obtained by the experimental measurements.

5. Static behavior

5.1. Discretization of the equations Eqs. (13)–(15) and relevant boundary conditions are discretized according to the Galerkin procedure. The first three eigenfunctions of the eigenproblem described in Section 4.3 are used as discretizing functions. The displacements vector {u1, u2, #3}T is expressed as a linear combination of the given z-function (eigenfunctions) and unknown t-functions qk(t):

8 9 2 9 38 0 0 /1 ðzÞ > > < q1 ðtÞ > < u1 > = = 6 7 /22 ðzÞ /23 ðzÞ 5 q2 ðtÞ u2 ¼ 4 0 > > > : : > ; ; 0 /32 ðzÞ /33 ðzÞ q3 ðtÞ #3

where the ith column of the matrix is the ith eigenfunction of the system. Eq. (21) are coupled in the u2 and #3 displacement components, since they reflect the coupled structure of the linear homogeneous equations of motion of the eigenproblem (16). The first three eigenfunctions of the beam are shown in Fig. 6, which resembles the structure of Eq. (21). A normalization rule that makes equal to one the eigenfunctions /22 and /23 at the free-end section is used. As a consequence the eigenfunctions /32 and /33 at the free-end section assume an absolute value much greater than one. The non-uniformity of the values assumed by the eigenfunctions is related to the fact that the functions /32 and /33 represent the torsional rotation of the beam while the functions /22 and /23 the transversal displacement of the beam. Finally, to obtain the discretized static equations, special care has been devoted to the torsional moment Mt applied at the free-end section. By referring to Fig. 4c, it reads: 0

Mt ¼ F  h ¼ F  h cos #3 ðlÞ

The experiment has been conducted by applying at the free-end section of the cantilever beam a torsional moment and measuring the displacements of some points of the two instrumented sections (Figs. 4 and 5). Comparison among the numerical results obtained by the different FE models and by the Galerkin discretization of the static equations, with and without the new nonlinear terms, is performed.

ð21Þ

ð22Þ

where the distance h0 between the forces F has been related to the deformed configuration of the free-end section of the beam. The effects associates with the loss of the parallelism of the forces F when the section rotates have been neglected since the distance between the loaded section and the pulleys is sufficiently great (see Fig. 5a and b).

Fig. 6. First three eigenfunctions and modal displacements of the free-end section.

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

1901

Fig. 7. Displacements measured on the Section 1: (a) displacement of point A; (b) displacement of point B; (c) displacement of point C with the new terms; (d) displacement of point C without the new terms (solid thick line: analytical model; circles: experimental test; triangles: FE models).

Fig. 8. Displacements measured on the Section 2: (a) displacement of point A; (b) displacement of point B (solid thick line: analytical model; circles: experimental test; triangles: FE models).

5.2. Comparison among numerical and experimental results A comparison among the displacement evaluated by the different models on the instrumented section is performed. From the results obtained by the discretized equations, by using Eq. (21) it is possible to obtain the displacements components of the shear center of the section u1, u2, #3. Then, making use of Eqs. (4) and (9), through Eq. (1) the displacement of a generic point on a generic section of the beam is evaluated. In the following, results obtained from the Galerkin discretization of the analytical equations will be labeled with ‘analytical model’ for the sake of brevity. These graphic conventions are always used: solid thick line refers to the results obtained by the analytical model; lines labeled with circles refer to the results obtained by the experimental test; lines labeled with triangles refer to the results furnished by the FEM models. In Fig. 7 results obtained by analytical, numerical (FEM2 – AnsysÒ) and experimental models are shown. In particular the displacement of the point A (Fig. 7a) and of the point B (Fig. 7b),

evaluated on the instrumented Section 1, are plotted versus the applied force F. It is possible to observe the good agreement among the models. Results obtained by the analytical model show a less rigid behavior than the other two. This fact can be explained reminding that in this model have been used flexural curvatures up-to quadratic terms (Eq. (9)1,2). If cubic terms are considered, even in the flexural curvature (Eq. (9)3), a general hardening of the beam could be observed as in Luongo et al. (1989) and the results would approach those obtained by FEM and experimental models. Very interesting is the comparison among results obtained taking and not taking into account the new nonlinear contributions. In Figs. 7c and d the displacement of the point C on the Section 1 with and without the nonlinear terms due to warping and torsional elongation are respectively shown. When the new term (underlined terms in Eqs. (13)–(15)) are considered, a good agreement among results obtained by the analytical, bi-dimensional numerical model (FEM2 – AnsysÒ) and experimental model is observed (Fig. 7c). On the contrary, when no new terms are taking

1902

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

into account, one-dimensional numerical model (FEM1 – AnsysÒ) and analytical model show a comparable behavior to each other, but very different to that furnished by the experimental model (Fig. 7d). In particular, in this last case, numerical and analytical models furnish displacements of point C opposite in sign with respect to the real case. In both cases, with and without new terms, the displacements of the point A and B evaluated with the different models are very close to each other. Finally, in Fig. 7c the dashed curve refers to the results obtained considering only the torsional elongation effects (neglecting nonlinear warping). As it is possible to observe, this curve is very close to the curve containing all the new contributions. This means that the main nonlinear effects are related to the torsional elongation terms, while only a small correction is related to nonlinear warping terms. The displacement of the point A and B evaluated on the instrumented Section 2 confirm what previously observed (Fig. 8). 6. Static bifurcation This analysis consists of the study of the flexural–torsional instability of the beam and on the comparison among results obtained by numerical models (FEM), analytical model and results of literature where critical values of the forces generating a flexural–torsional instability are already evaluated for classical opencross section beams. 6.1. Discretization of the equations In this case the equations and the relevant boundary conditions are discretized by using different functions with respect to the previous case. The eigenfunctions that are uncoupled in the displacements, given by Eq. (16) when x1C ¼ 0, corresponding to a double symmetric open cross-section thin-walled beam, are used. This choice is related to the fact that these uncoupled eigenfunctions seem to be able to better analyze the uncoupled flexural solutions and their stability. This is confirmed by a convergence analysis where it has been shown that a smaller number of eigenfunctions is required to better describe the stability of the two pure flexural deflections occurring when the free-end section is loaded by a force acting on the shear center and directed along the x1 or the x2 axes. The independent displacements vector {u1, u2, #3}T is then expressed:

9 8 9 2 38 /11 ðzÞ 0 0 > > < q11 ðtÞ > = < u1 > = 6 7 /21 ðzÞ 0 u2 ¼ 4 0 5 q21 ðtÞ > > > : ; : > ; 0 0 /31 ðzÞ q31 ðtÞ #3 2 6 þ4

first three eigenfunctions

/12 ðzÞ 0 0

9 8 9 38 > < q12 ðtÞ > = > < /1s > = 7 /22 ðzÞ 0 5 q22 ðtÞ þ /2s > > > : ; > : ; 0 /32 ðzÞ q32 ðtÞ 0 0

0

ð23Þ

second three eigenfunctions

where /ij(z) (j = 1, 2) are the first two uncoupled eigenfunctions for each component and /is are static functions used to represent the

behavior of a concentrated force as experienced for different problems in Biondi and Muscolino (2005) and Luongo and D’Annibale (submitted for publication). The additions of /is make it possible to reduce the number of the discretizing functions necessary to obtain the correct bifurcation points and the post-critical behavior. The positions /2s ¼ 0 or /1s ¼ 0 must be adopted when the freeend section is loaded by a force applied on the shear center and respectively directed along the x1 or the x2 axes. The static functions read:

^i  /is ¼ /

2 X

^ik /ik q

ði ¼ 1; 2Þ

ð24Þ

k¼1

^ i are the static deflection due to the forces applied on the where / shear center of the free-end section and directed along the xi axis; P ^ k /ik qik represents the projection of this deflection into the base of the used eigenfunctions (i = 1, 2). By referring to a force acting ^ 1 and the static funcalong the x1 axis, in Fig. 9 the static deflection / tion /1s are shown. Even if the static function /1s is very small, it is able to describe the presence of the concentrated load at the free end better than the use of only eigenfunctions (see the evolution of the shear function EI2 /000 1s in Fig. 9). 6.2. Numerical results First the stability of the pure flexural deflection along the symmetry axis x1 is discussed (Fig. 2a). The increasing force F applied on the shear center of the free-end section of the beam, directed along the x1 axis, plays the role of bifurcation parameter. Critical and post-critical behaviors of the beam are studied by solving the discretized nonlinear static equations of the beam by using AutoÒ (Doedel and Oldeman, 2009). In Fig. 10 bifurcation diagrams are shown, where stable paths are indicated with thick solid lines, while unstable paths with thin solid lines. The post-critical paths of the six modal displacements qij (Eq. (23)) are reported. At point A the first static bifurcation manifests itself (F = 1210 N). The modal components of the displacement u2 and #3 (Figs. 10b and c), that are zero before the bifurcation, show a classical post-critical stable behavior. The symmetry of the bifurcated paths is strictly related to the symmetry of the section with respect to the x1 axis. In Figs. 10d and e the two possible symmetric post-critical configurations of the free-end section are reported. Finally, with dashed lines, perturbed paths are shown. They refer to the cases where the force F is not exactly applied on the shear center of the section, but it has an eccentricity equal to 2% of h (±0.001 m). As expected, they are very close to the post-critical paths and confirm the critical behavior. The stability of the pure flexural deflection along the axis x2 is also studied (Fig. 2a). The increasing force F applied on the shear center of the free-end section of the beam, directed along the x2 axis, plays the role of bifurcation parameter. In Fig. 11 bifurcation diagrams are shown and the post-critical paths of the six modal displacements qij are reported. At points A and B two static bifurcations manifest themselves (FA = 1700 N, FB = 2104 N).The modal components of the displacement u1 and #3 (Fig. 11b and c), that are zero before the bifurcation, show a non-symmetric post-critical

Fig. 9. Static function in the x1 direction.

1903

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

(a)

(b)

(c)

(d)

(e)

Fig. 10. Bifurcation diagrams: (a) modal component of the displacement u1; (b) modal component of the displacement u2; (c) modal component of the rotation #3; (d, e) the two symmetric flexural–torsional deflection after bifurcation point A.

Fig. 11. Bifurcation diagrams: (a) modal components of the displacement u2; (b) modal components of the displacement u1; (c) modal components of the rotation #3; (d) flexural–torsional deflection after bifurcation point A; (e) flexural–torsional deflection after bifurcation point B.

behavior due to the non-symmetry of the section with respect to the x2 axis. In Figs. 11d and e the two possible post-critical configurations of the free-end section are reported. They are not symmetric and can occur after the two critical loads FA and FB,

respectively. In particular the configuration described in Fig. 11e after the bifurcation point B can be never reached, since the postcritical behavior after the point A is unstable. Finally, with dashed lines, perturbed paths are shown. They refer to the cases where the

1904

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

Table 4 Critical loads. One-dimensional ThreeAnalytical Literature (N) model model (FEM1) (N) dimensional model (FEM2) (N) (AM) (N) ADINAÒ ANSYSÒ ADINAÒ ANSYSÒ 2620

2570

1350

1290

1210

2568

1020

1020

1690

1650

1700

1017

force F is not exactly applied on the shear center of the section, but it has an eccentricity equal to 1% of h (±0.0005 m). Also in this case, they are very close to the post-critical paths previously discussed. In Table 4 the values of the first critical loads, obtained by different models, are reported. It is interesting to observe that the FE models based on a one-dimensional thin-walled beam element (FEM1) are in good agreement with the results of the scientific and technical literature (Pignataro et al., 1990). However, these values are very different from the results obtained by the more refined FE models based on bi-dimensional shell finite elements (FEM2). The fact that the results obtained by the analytical model (AM) are very close to those obtained by the FEM2, confirms the correctness of the proposed formulation and the importance of the new nonlinear contributions, previously described. In the classical nonlinear open-cross section beam models the flexural–torsional instability is strictly connected to the nonlinear terms of the equations generated by the cubic terms of the flexural curvatures. In this model it is proven that the new nonlinear terms related mainly on torsional elongation are prevailing with respect to the previous cited terms of the classical beam models, since they appear in the equations at a lower perturbation level. Also in the stability analysis the absence of the torsional elongation does not permit to correctly evaluate the critical load. When the nonlinear warping is considered, the correction of the critical loads obtained by taking into account only the torsional elongation is small. This circumstance explains why this model furnishes critical loads shown in the AM column of Table 4 that drastically change when the new terms are not taken into account (FEM1 and Literature columns of Table 4). Finally, some consideration must be given to justify the discretization adopted for the bifurcation analysis. For the study of the flexural–torsional instability of the pure flexural solution along the x1 axis the first three uncoupled eigenfunctions (without static functions /1s ) could be sufficient. On the contrary the instability of the pure flexural solution along the x2 axis requires a higher number of discretizing functions. The first six uncoupled eigenfunctions together with the static function /2s are the minimum number of discretizing functions able to better evaluate the critical and the post-critical behavior of the beam. If no static function is used, at least nine uncoupled eigenfunctions are necessary to achieve the goal. If instead coupled eigenfunctions are used, more than six functions need to obtain acceptable results. For uniformity, in both the bifurcation analyses performed, six uncoupled eigenfunctions plus a static function have been always used.

developed (Di Egidio et al., 2003, Part I), has been specialized to describe the behavior of a cantilever beam with a monosymmetric section when the torsional curvature is greater than the flexural ones. The simplified static equations obtained make it possible to stress the role of the nonlinear terms due to nonlinear warping and torsional elongation of the longitudinal fibers. A preliminary study has been developed to determine the geometrical characteristic of the channel cross-section of technical interest for which the ratio between torsional and flexural curvatures is high. The first analysis performed concerns the behavior of the beam under static loads. It has been conducted in a laboratory on an experimental model and by means of nonlinear finite element models. In the second analysis the stability of the static equilibrium branches has been studied by using only different numerical models (FEM). Since the first scope of the paper is the validation of the analytical model, a comparison among the results of this one and those obtained by numerical (FEM) and experimental models has been carried out. The comparison has shown the important role of the nonlinear warping and the effectiveness of the analytical model both in the static and in the stability analyses and the important role played by nonlinear warping and torsional elongation. In particular the model furnishes values of critical loads in flexural–torsional stability that the classical nonlinear one-dimensional beam models are not able to describe correctly. From a technical point of view it has been shown that the literature results for the critical loads are very far from the results obtained by means of the refined model proposed. In one case – flexural–torsional buckling along symmetry cross section axis – is almost double the correct one, i.e. not in favor of safety, and in another case – flexural–torsional buckling along no symmetry axis – is almost half the correct one, in favor of safety. Acknowledgements The authors recognize the contribution of the Prof. A. Luongo in the developing of the analytical model (Di Egidio et al., 2003, Part I; Vestroni et al., 2006). This work was partially supported under the FY 2007-2008 PRIN Grant. Appendix A The warping functions in Eq. (5) are defined as:

~1 ¼ b  b 1  k12 x2 b 1  ~ 2 þ k12 x1 b2 ¼ b2  b ~3 ¼ b  b 3  k31 x1  k32 x2  2b ~6 k45 b

ðA:1Þ

3

1 ; k12 ; k31 ; k32 are: where b

R k3i ¼

A

b3 xi dA ; I2

R k12 ¼

A

ðb2 x1  b1 x2 ÞdA ; ðI1 þ I2 Þ

R j ¼ b

A

bj dA A

ðA:2Þ

and where bj read:

Z 1 c 2 r ðcÞcoswðcÞdc 2 0 Z c 1 r 2 ðcÞsinwðcÞdc b2 ðcÞ ¼  2 0 Z c b3 ðcÞ ¼ ½a1 ðcÞrðcÞ þ 2b1 ðcÞ cos wðcÞ þ 2b2 ðcÞsinwðcÞdc

b1 ðcÞ ¼ 

ðA:3Þ

0

7. Conclusion

Appendix B

A nonlinear one-dimensional model of a thin-walled, open cross-section, shear and axially-undeformable beam, previously

The geometrical and the mechanical characteristics of the beam are reported in Tables 1 and 2.

A. Di Egidio, F. Vestroni / International Journal of Solids and Structures 48 (2011) 1894–1905

References Ascione, L., Feo, L., 1995. On the mechanical behaviour of thin-walled beams of open cross-section. Int. J. Eng. Anal. Des. 2, 14–34. Biondi, B., Muscolino, G., 2005. New improved series expansion for solving the moving oscillator problem. J. Sound Vib. 281, 99–117. Brincker, R., Zhang, L., Andersen, P., 2001. Modal identification of output-only systems using frequency domain decomposition. Smart Mater. Struct. 10, 441– 445. Crespo da Silva, M.R.M., 1988. Non-linear flexural–flexural–torsional–extensional dynamics of beams – I. Formulation. Int. J. Solids Struct 24 (12), 1225–1234. Crespo da Silva, M.R.M., 1991. Equations for nonlinear analysis of 3D motions of beams. Appl. Mech. Rev. 44 (11), 51–59. part 2. Crespo da Silva, M.R.M., Glynn, C.C., 1978. Nonlinear flexural–flexural–torsional dynamics of inextensional beams – I. Equations of motion. J. Struct. Mech. 6 (4), 437–448. Crespo da Silva, M.R.M., Zaretzky, C.L., 1994. Nonlinear flexural–flexural–torsional interactions in beams including the effect of torsional dynamics – I. Primary resonance. Nonlinear Dyn. 5, 3–23. Doedel, E.J., Oldeman B.E., 2009. AUTO-07: continuation and bifurcation software for ordinary differential equations. Technical Report, Department of Computer Science, Concordia University, Montreal, Canada. Available from . Di Egidio, A., Luongo, A., Vestroni, F., 2003a. A nonlinear model for open crosssection thin-walled beams – Part I: Formulation. Int. J. Non-Linear Mech. 38 (7), 1067–1081. Di Egidio, A., Luongo, A., Vestroni, F., 2003b. A nonlinear model for open crosssection thin-walled beams – Part II: Forced motion. Int. J. Non-Linear Mech. 38 (7), 1083–1094.

1905

Ghobarah, A.A., Tso, W.K., 1971. A non-linear thin-walled beam theory. Int. J. Mech. Sci. 13, 1025–1033. Luongo, A., Rega, G., Vestroni, F., 1989. Non resonant non-planar free motions of inextensional non-compact beams. J. Sound Vib. 134 (1), 73–86. Luongo, A., D’Annibale, F., submitted for publication. Linear stability analysis of multiparameter dynamical systems via a numerical-perturbation approach. AIAA J. Mollmann, H., 1982. Finite displacements of thin-walled beams. Part 1, geometrical foundations. Report No. 252, Danish Center for Appl. Math. and Mech., Technical University of Denmark, pp. 1–48. Pignataro, M., Rizzi, N., Luongo, A., 1990. Bifurcation, Stability and Postcritical Behaviour of Elastic Structures. Elsevier Science Publishers, Amsterdam. Rosen, A., Friedmann, P., 1979. The nonlinear behavior of elastic slender straight beams undergoing small strain and moderate rotations. J. Appl. Mech. 46, 161– 168. Sapountzakis, E.J., Mokos, V.G., 2003. Warping shear stresses in nonuniform torsion by BEM. Comput. Mech. 30, 131–142. Sapountzakis, E.J., Tsiatas, G.C., 2007. Flexural–torsional vibrations of beams by BEM. Comput. Mech. 39 (4), 409–417. Sapountzakis, E.J., Tsiatas, G.C., 2008. Nonlinear inelastic uniform torsion of bars by BEM. Comput. Mech. 42 (1), 77–94. Schulz, M., Filippou, C., 1998. Generalized warping torsion formulation. J. Eng. Mech. 124, 339–347. Vestroni, F., Di Egidio, A., Luongo, A., 2006. Nonlinear warping and torsional elongation effects on the response of an open cross section beam. In: Proceedings of the 6th European Solid Mechanics Conference – ESMC 2006, Budapest, Hungary, 28 August–1 September. Vlasov, V.Z., 1961. Thin-walled Elastic Beams. Monson, Jerusalem.