Elasticity-based free vibration of anisotropic thin-walled beams

Thin-Walled Structures 95 (2015) 73–87

Contents lists available at ScienceDirect

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

Elasticity-based free vibration of anisotropic thin-walled beams Paul R. Heyliger Department of Civil and Environmental Engineering, Colorado State University, Fort Collins, CO 80523, United States

art ic l e i nf o Article history: Received 24 December 2014 Accepted 19 June 2015

a b s t r a c t The free-vibration of anisotropic open and closed thin-walled beams is studied using approximate solutions to the three-dimensional equations of motion of a linear elastic solid. The three displacement components are approximated using grouped polynomials over both the cross-section coordinates and independent polynomial or trigonometric functions over the axial coordinate of several beams with representative boundary conditions. This method is applied to several sections with varying degrees of anisotropy and shape with both isotropic and orthotropic material properties. Results are compared with those computing using beam or shell theories, and several new results are presented for these geometries. & 2015 Elsevier Ltd. All rights reserved.

1. Introduction Three-dimensional solids that have one-dimension much larger than the other two and a cross-section composed of elements of that are generally long and thin are usually classified as thinwalled beams. Their behavior, which can combine elements of motion usually associated with axial, flexure, and torsional motion, has been extensively studied ever since the pioneering work of Wagner and Pretschner [1], Bleich and Bleich [2], and Kappus [3]. As noted by Gjelsvik [4], the comprehensive theory of open thinwalled sections of Vlasov was also developed during the same period as these early works but not translated into English until 1963 [5]. By this time Goodier [6] and Timoshenko [7] had independently developed general theories. Additional work by Gere and Lin [8] on coupled vibration behavior outlined some of the more significant issues related to these components. Numerous computational models have been further developed, including the dynamic stiffness matrix method of Banerjee and co-workers [9] and the finite element model of Mei [10]. More recent applications to composite beams have also been considered by Kollar [11], Cortinez and Piovan [12], Piovan and Cortinez [13]. Vibration of thin-walled beams under initial stress has been studied by Machado and Cortinez [14]. Benscoter [15] developed an early theory for thin-walled beams with closed cross-section, where the general behavior becomes more complicated. Significant progress has been made since this time, including a number of approaches for anisotropic beams, including those of Hodges and co-workers [16] and Song and Librescu [17]. Summaries of much of this work have been given by E-mail address: [email protected] http://dx.doi.org/10.1016/j.tws.2015.06.014 0263-8231/& 2015 Elsevier Ltd. All rights reserved.

Murray [18], Librescu and Song [19], and Hodges [20] for beams of open and closed section under both static and dynamic response. By far the most comprehensive studies of thin-walled beams have used beam, plate, or shell theories in which displacement fields are assumed over the domain of the plate that restrict the motion of the displacement components [21]. By their own definition, the displacement fields of all beam theories, thin-walled or not, require some sort of kinematic restriction. This usually takes the form of assumed components of representative displacement and/or stress over the cross-section of the beam. This is intended to allow a rapid analysis that captures the key mechanics of these solids. Approximate beam models are usually compared with either more detailed elasticity results or experimental evidence to assess their accuracy. It is less typical to allow for a full representation of the usual displacement components consistent with the three-dimensional theory of elasticity. In this study, elasticity-based displacement fields are assumed over the crosssection of both isotropic and anisotropic beams. This is accomplished by applying the Ritz method directly to the three-dimensional equations of elasticity using power-series approximations in the form introduced by Visscher [22]. The intent of this work is to determine the level of accuracy of more simplified models, provide an alternative means of analysis for these solids, and present numerical results for thin-walled beams to provide a basis for future comparison.

2. Geometry and governing equations We consider a three-dimensional solid whose cross-section coordinates are defined in the (x1,x3) or (x–z) plane with a much

74

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

larger dimension in the x2 or y-axis with a length of L. The beam is assumed to be composed of an anisotropic material whose principal material axes are aligned with the (x,y,z) axes. This material is assumed to be originally orthotropic and having gone through a rotation in the x–y plane at an angle θ as measured from the positive x-axis (the horizontal component of the cross-section axis). The beam is hollow with a uniform wall thickness of t. 2.1. The weak form The governing equations used for this study are the three-dimensional equations of linear elasticity with an anisotropic constitutive tensor that can represent what amounts to a monoclinic material. These are not solved explicitly at each point in the domain, but instead approximate solutions are sought for their weak form as expressed within Hamilton's Principle [23]. This is usually expressed as t

0=− +

∫0 ∫V {σ1δ ϵ1 + σ2δ ϵ2 + σ3δ ϵ3 + σ4δ ϵ4 + σ5δ ϵ5 + σ6δ ϵ6} dV dt 1 δ 2

t

∫0 ∫V ρ(u̇2 + v2̇ + ẇ 2) dV dt

(1)

Here the conventional contracted notation has been used for Cauchy stress components (s11 ¼s1, s23 ¼s4, ϵ1 = ϵ11, ϵ4 = 2ϵ23, C1111 ¼C11, C1123 ¼C14, etc.) and it is understood that the 1, 2, and 3 directions are (x1 ¼x, x2 ¼y, and x3 ¼z). The general constitutive law that links the Cauchy stress to the components of infinitesimal strain ϵij can be written as

σij = Cijkl ϵkl

(2)

where the strain components are linked with the three displacement components u, v, and w in the x, y, and z directions as

∂u ∂v ∂w , ϵ2 = ϵ22 = , ϵ3 = ϵ33 = ∂x ∂y ∂z ∂v ∂w ∂u ∂w ∂v ∂u , ϵ5 = γ13 = , ϵ6 = γ12 = . ϵ4 = γ23 = + + + ∂z ∂y ∂z ∂x ∂x ∂y (3) ϵ1 = ϵ11 =

The general constitutive relation for the materials used in this study can be expressed in matrix form as

⎧ ⎫ ⎡ C11 ⎪ σ1⎪ ⎢ ⎪ σ ⎪ ⎢C12 2 ⎪ ⎪ σ3 ⎪ ⎪ ⎢C ⎨ ⎬ = ⎢ 13 σ ⎪ 4⎪ ⎢ 0 ⎪ σ5 ⎪ ⎢ ⎪ σ6 ⎪ ⎢ 0 ⎪ ⎩ ⎪ ⎭ ⎢⎣C16

C12 C13

0

0

C22 C23

0

0

C23 C33

0

0

0

0

C44 C45

0

0

C45 C55

C26 C36

0

0

C16 ⎤⎧ ⎫ ⎥⎪ ϵ1⎪ C26 ⎥⎪ ϵ ⎪ ⎪ 2⎪ C36 ⎥⎥⎪ ϵ3 ⎪ ⎨ ⎬ 0 ⎥⎪ ϵ4 ⎪ ⎥⎪ ϵ ⎪ 0 ⎥⎪ 5 ⎪ ϵ6 C66 ⎥⎦⎪ ⎩ ⎪ ⎭

(4)

These stiffness components are associated with either a naturally occurring monoclinic material with 13 independent elastic constants or a rotated orthotropic material. In the latter case, the final Cij as listed in the 6  6 matrix are a function of the on-axis elastic constants and the angle of rotation of originally orthotropic properties about the 3 or z-axis. Their explicit forms are given in Reddy [24]. 2.2. The Ritz model The computational approach used in this study relaxes many of the assumptions used in most beam theories. Rather than introducing specific displacement fields over the beam domain, the three point-wise displacement components are approximated directly using the three-dimensional equations of elasticity. There is no need to represent any rotational variables or their equivalent in this theory, since each of the displacements is computed at every location in the beam.

The general form for the displacements can be expressed as N

ui (xj , t ) = ϕoui(xj ) +

∑ Aip(t )ϕpui(xj) p=1

(5)

Here the Aip represent unknown constants that depend on time. The three values of i indicate that there are three components of displacement and that three independent sets of functions can be used to approximate each of those displacements. The ϕo terms represent the simplest functions that satisfy the essential boundary conditions for that displacement direction. In the study discussed in this paper, the initial boundary conditions will be assumed to be zero, with no initial displacement or velocity. The ϕo terms are therefore all equal to zero for the considered case. In these approximation equations, ϕp represents a selected approximate function for each respective direction. These functions must satisfy general requirements of independence as outlined by Reddy [23]. By using a large number of these approximation terms for each displacement component, very accurate solutions can be determined for the mode shapes and frequencies of vibration for a beam. In this study, polynomial series are used over the cross-section of the beam and then combined with other polynomials or transcendental functions over the axial length of the beam. This class of approximation has a long and rich history in the field of mechanics. Among the first to use it for the case of vibrating solids was Demarest [25], who used Legendre polynomials over the three Cartesian coordinate directions to find the natural modes of vibrating parallelepipeds. This method was further studied by Ohno [26], who also discussed the splitting of the resulting eigenvalue problem into a smaller number of equivalent eigenvalue problems using group theory. Polynomial expansions were also used in general cylinder [27] and sphere [28] vibrations. A significant and powerful expansion of this approach was introduced by Visscher [22], who expanded the three displacement components in terms of power series regardless of the shape of the solid. This is generally written as

ϕpu(x, y , z ) = xi y j z k

(6)

This is a very attractive family of approximations in that it frequently yields very simple integration results over the domain of the solid. The nature of this integration will be discussed for specific cases in the sections that follow. This type of approximation has seen extensive application in resonating solids [29], vibrating beams [30], and related areas such as granular media [31]. In terms of applications to solids with more general beam geometry, Carrera and co-workers have recently developed a unified beam theory using a modification of this type of approximation over the beam cross section and have extensively applied it to a wide array of problems [32,33]. If the displacement constants cip are replaced with the letters c, d, and e for the three displacements (1 ¼u, 2¼ v, 3 ¼w) for any index p, it is possible to generalize the matrix form of the Ritz model applied to Hamilton's principle following the assumption of general harmonic motion. Specifically, the Ritz coefficients can be written as

A1p (t ) = cpsinωt

(7)

A2p (t ) = dpsinωt A3p (t ) = epsinωt

(8)

Substitution of these conditions and the general Ritz approximation into Hamilton's principle leads to the generalized eigenvalue problem that can be represented in matrix form as

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

⎡[M11] 0 ⎡[K11] [K12] [K13]⎤⎧ ⎫ ⎫ 0 ⎤⎥⎧ {c} {c} ⎪ ⎢ ⎥⎪ ⎢ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 2 22 21 22 23 ⎢[K ] [K ] [K ]⎥⎨{d}⎬ = ω ⎢ 0 [M ] 0 ⎥⎨{d}⎬ ⎥⎪ ⎪ ⎢ ⎥⎪ ⎪ ⎢ 0 [M33]⎦⎪ ⎣ 0 ⎣[K 31] [K23] [K 33]⎦⎪ ⎭ ⎩ {e} ⎪ ⎭ ⎩ {e} ⎪

(9)

The elements of the stiffness ([K ]) and mass ([M ]) matrices are given in the Appendix. In the specific cases that follow, the (x,z) plane is taken as the cross-section of the beam and the y-axis is the longitudinal axis of the beam. For the most part, the approximation functions in the ydirection appear as power series. For certain types of support conditions this may require some adjustment or replacement. This will be discussed for beams under simple support studied in the next section. Computationally, using global power series over the crosssection greatly simplifies the integration required for the relevant coefficient matrices. For hollow sections, for example, the integration can be completed first over the outer bounding dimensions of the solid and then subtracting the same integral evaluated over the inner bounding surfaces of the cross-section. The operations for rectangular sections are straightforward. For the elliptical section, the integration requires a trigonometric substitution that can be evaluated in closed form for any powers of x and z [31]. One additional advantage of using global approximation functions over the beam cross-section and the length is the potential to group the specific approximations according to the physical character of the displacement components they are being used to represent [26,34]. Various conditions of cross-sectional symmetry, axial wavelength, and grouping of the various approximation functions can drastically reduce the size of the matrix eigenvalue problem. Rather than generalizing this discussion, these features are discussed for the specific geometries for which they are used in the examples that follow. 2.3. Beams with simple support: semi-analytic FEM model In cases where the axial (or y) dependence of the displacement functions are known, the form of the approximation functions can be modified to reduce the dimension of the problem and fix the axial dependence as a given quantity. In other words, when the axial dependence is fixed the problem becomes two-dimensional in the cross-sectional variables (x,z). This allows for the introduction of a second independent method of analysis with which results of the present model can be compared. In the case of some beams under the end conditions of simple support beam, most beam theories usually require that the components of transverse displacement are zero at the end-points and also that the resultant moment quantities are equal to zero. Using the elasticity theory of the present model, identical displacement conditions are imposed at the endpoints of the beam for all points over the cross-section at all times. This implies that u(x, 0, z, t ) = u(x, L, z, t ) = w (x, y, 0, t )= w (x, L, z, t ) = 0 . In addition, the traction vector component along the axis ty ¼syy ny ¼0. For the types of materials considered in this study, each of these conditions is satisfied for beams up to and including orthotropic symmetry (C16 ¼C26 ¼ C36 ¼C45 ¼0) using the approximations n

u(x, y , z, t ) =

∑ Uj(t )Ψ ju(x, y)sinky j=1 n

v(x, y , z, t ) =

∑ Vj(t )Ψ jv(x, y)cosky j=1 n

w (x, y , z, t ) =

∑ Wj(t )Ψjw(x, y)sinky j=1

(10)

Here k is the wave number and the capital letters for the field

75

variables are each spatial constants associate with the individual subscript j, each of which corresponds to independent shape functions as denoted by Ψ and each of which can still vary with time. When combined with conventional finite element approximations in two dimensions, they represent the nodal values of the displacement components. As for the general Ritz method, the superscript on each function indicates that the approximation functions for each of the variables need not be the same. Depending on the geometry of the problem each approximation function can take a specific form. The variations of the primary field variables within Hamilton's Principle are then taken as

δu = Ψiu(x, y),

δv = Ψiv(x, y),

δw = Ψi w(x, y)

(11)

Although solutions can be determined for the general transient case, the focus in this study is on the calculation of periodic frequency for a given wave number k. Hence it is further assumed that the time dependent behavior for the spatial constants is in the form given by general periodic motion as

Uj(t ) = Ujsinωt ,

Vj(t ) = Vjsinωt ,

Wj(t ) = Wjsinωt

(12)

Here t is the time and ω is the natural frequency of free vibration. Substitution of these functions and collecting terms allows the weak form from Hamilton's principal to be expressed in matrix form in a manner nearly identical to the original Ritz method:

⎡[M11] 0 ⎡[K11] [K12] [K13]⎤⎧ ⎫ ⎫ 0 ⎤⎥⎧ {U} ⎪ {U} ⎪ ⎪ ⎢ ⎥⎪ ⎢ ⎪ ⎪ ⎪ ⎪ 2 22 21 22 23 ⎢[K ] [K ] [K ]⎥⎨ {V} ⎬ = ω ⎢ 0 [M ] 0 ⎥⎨ {V} ⎬ ⎪ ⎪ ⎥⎪ ⎢ ⎥⎪ ⎢ 0 [M33]⎦⎪ ⎣ 0 ⎣[K 31] [K23] [K 33]⎦⎪ ⎩{W}⎪ ⎭ ⎩{W}⎪ ⎭

(13)

There are two primary differences between this form and that of the general Ritz model that is the focus of the present study: (1) the integration needs to occur only over the beam cross-section and (2) the constants (in this case the nodal variables) exist only over a single cross-section. This drastically reduces the computational size of the problem. This model is used as a method of comparison for the beams with simple support in the sections that follow.

3. Results and discussion In the results that follow, the present method is labeled 3DRE to represent the three-dimensional Ritz elasticity structure of the analysis. This is used in all tabulated results for the examples considered. Several general shapes are considered. The first are somewhat standard beam cross-sections that are not necessarily thin-walled but are useful to study the accuracy of the present 3DRE model. 3.1. General sections: the reduction to basic beam theories One feature of Ritz-based methods based on power-series approximations of the displacement components is that it allows a simple reduction to capture the behavior of many elementary theories. Han et al. [35] have published a complete assessment of the Euler–Bernoulli, Rayleigh, and Timoshenko beam theories under a variety of boundary conditions. These results can then be directly compared to the current elasticity-based Ritz models to show the agreement between approaches. In a power-series representation of the displacement vector, the displacement dependence on the cross-section coordinates x and z can be restricted to exactly match a number of beam theories for sections that are symmetric about either of these axes (it is assumed for this section only that the cross-sections are symmetric about the z-axis): The Rayleigh theory can be recovered by

76

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

setting the Poisson ratio used to calculate the elastic stiffnesses Cij to zero, setting i ¼k¼ 0 as a single term approximation for the displacement components u and w for all y, and keeping only the approximation i ¼0 and k¼ 1 as the approximation for v in Eq. (6). In addition, the shear modulii are increased by 5–6 orders of magnitude to reduce the influence of shear deformation. The Timoshenko theory can be simulated by using the same restrictions as for the Rayleigh theory but using the exact shear modulii to account for first-order shear deformation and using an appropriate shear correction coefficient. The clamped–free beam is a very common configuration that has been extensively studied for a variety of geometries. Han and co-workers [35] have given formulas for the characteristic frequency equations. The frequencies can be computed for any range of input parameters. For the Rayleigh theory, these include the modulus of elasticity E, the density ρ, the length L, and the crosssectional area A. For the Timoshenko model these parameters are required in addition to the shear modulus G and the value of the shear coefficient k. This latter value has been the subject of intense study but there are fairly well accepted values for typical crosssections and noted by Han [35]. As an example, a cantilever beam was studied whose relative dimensions were 1  1  10 in the cross-section and length directions. Two shapes were considered: the square and the circle. The Rayleigh and Timoshenko approximations were used over the cross-section and so the Poisson ratio was assumed to be zero. Along the axis of the beam, power series in y (y1, y2, and so on) were used for each of the displacement components to denote the requirement that all displacements are equal to zero at the clamped end where y¼ 0. For each cross-sectional shape, two additional refinements were considered: the solid section and the thin-walled closed section where the wall thickness is 1/100 of the outer dimension. The results gave agreement within four or five figures between the constrained elasticity theory developed here and the Rayleigh and Timoshenko theories. The values of the shear coefficient [36,37] were taken from Han [35] who had summarized the work of Cowper [38]. 3.2. Open sections The specific open thin-walled sections considered are the unsymmetric channel shown in Fig. 1 and the symmetric semicircle

Fig. 1. The beam cross-section of the doubly unsymmetric channel. For the results presented here, the origin of the (x,y) system was located at the lower-left interior corner of the beam cross-section. This required integration of the coefficient matrices over four rectangular subdomains.

of Fig. 2. Both of these sections have seen extensive study, and the goal of their inclusion here is to evaluate the accuracy of earlier theories and models within the context of the present elasticity solution. 3.2.1. The doubly unsymmetric channel The doubly symmetric section shown in Fig. 1 was originally studied by Yaman [39] and has been extensively used by a number of other researchers. It is an excellent example of a doubly symmetric section as there are no mirror planes that intersect the x–y cross-sectional plane. The material properties used were an elastic modulus E of 70 GPa, a shear modulus G of 26 GPa, and a density of 2700 km/m3. Yaman originally gave only the integrated section properties rather than the actual cross-sectional dimensions. The latter are needed to integrate the power series expansions using the present elasticity model. Since these properties were not given by Yaman, the integrated quantities were used to compute the leg dimensions and thickness by minimizing the absolute difference between the area and the elements of the inertia tensor. Using this approach, the Yaman beam was determined to have a constant wall thickness of t¼t2 ¼ 1.275 (all mm), d1 ¼13.29675, d3 ¼25.9616, and d2 ¼39.375. This gives A ¼9.68 (9.68)  10  5 m2, Ixx ¼5073.64 (5080)  10  6 m4, Iyy ¼22 400 (22 400)  10  6 m4, and Ixy ¼4256 (4250)  10  6 m4 where the values used by Yaman have been given in the parentheses. This gives a total absolute difference of 13.23  10  6 m4 for the 3 values of the inertia tensor [40]. The level of agreement between competing formulations with regard to this cross-section has been somewhat mixed. Ambrosini [41] has provided a recent summary and discussion regarding the results of several researchers who have considered this crosssection, including the work of Arpaci and Bozdag [42] and Tanaka and Bercin [43]. Much of this discussion centers around the use of various elements of the inertia tensor of the beam cross section and other features related to the beam theories used to study the beam behavior. This discussion does not apply to the present method since the actual cross-section is used along the three-dimensional elasticity tensor components. But Ambrosini [41] notes the differences in the lower frequencies among these various studies. To assess the overall accuracy of the present method, the lowest four frequencies of the simply supported beam were computed and compared with the values of Ambrosini [41] and those of Yaman [39] along with the results from the two-dimensional semi-analytic finite element model and a fully three-dimensional finite element model. In the case of the two-dimensional finite element model, the wave number k can be directly input to extract the appropriate natural frequencies under the conditions of simple support. For the lowest four modes, three of the frequencies are associated with the wave number k ¼ π. This is consistent with modal vibration patterns that possess variations

Fig. 2. The beam cross-section of the cantilever semi-circle. The origin of the (x,y) system is located at the center of the full-circle extension of this geometry.

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

along the axial length that are symmetric about the midpoint of the beam. The third lowest frequency is associated with the lowest anti-symmetric mode where k ¼2π. By anti-symmetric we mean about the mid-point of the beam. The mesh used to represent the cross-section for both the two and three-dimensional finite element models has, along the three legs, 4, 12, and 8 element divisions and four element divisions through the thickness. The two corner regions therefore have 4  4 meshes. A total of 128 elements with 165 nodes are used to model the complete crosssection. The nature of the semi-analytic approximation ensures that all boundary conditions are explicitly satisfied at the endpoints of the beam. For the three-dimensional model, this same cross-section is used along with 40 element divisions along the length of the beam. Symmetry conditions were used so that only one-half of the beam length is modeled. Hence the boundary conditions at the simply supported end and the mid-point of the beam for the symmetric modes of the three-dimensional model are given by

u(x, y , 0) = v(x, y , 0) = σzz (x, y , 0) = 0

(14)

w (x, y , L/2) = σxz (x, y , L/2) = σyz (x, y , L/2) = 0

(15)

As is usually the case in standard finite element formulations, the boundary conditions on the stress components are not satisfied pointwise but rather in an integral sense over the entire crosssection of the beam. The results of these analyses are given in Table 1. There is good agreement between the various methods for the lowest frequency and also good agreement between the present methods, but significant variations in the remaining frequencies. The modal shapes for the half-beams are also shown in Fig. 3a–d. 3.2.2. The thin-walled semi-circle Friberg [44] was the first to study the vibration of the cantilevered beam shown in Fig. 2 using Vlasov beam theory. Since that time numerous researchers have also studied this particular crosssection, including Jun [45,46], de Borbon and Ambrosini [47], and Petrolo and co-workers [48]. This beam section has some advantages over the previous section in that the vibrational modes for isotropic beams can be split into groups that reflect symmetry or asymmetry about the vertical axis that bisects the section. These two modal groups have displacement components that are either even or odd about this axis. For symmetric modes, u is odd while v and w are even. For antisymmetric modes, u is even while v and w are odd about x. The properties of this section include an inner diameter di ¼45 mm, a wall thickness of 4 mm, a length of 820 mm, an elastic modulus of E ¼68.9 GPa, a shear modulus of 26.5 GPa, and a density of 2711 kg/m3. Two types of support conditions were considered: the fixed– free beam originally studied by Friberg [44] and the simply supported beam of the same cross-section considered by Jun [45,46] and de Borbon and Ambrosini [47]. In the present model, 3 terms were used in the x and y directions for the three variables (these terms include only those represented by the appropriate Table 1 Frequencies of the doubly unsymmetric channel. Mode

1 2 3 4

Frequency 3D FEM

2D FEM

[41]

[39]

3DRE

50.59 66.98

47.38 (k ¼ π) 64.60 (k ¼ π) 144.5 (k ¼ 2π) 158.5 (k ¼π)

51.39 96.48 188.84 204.81

51.87 114.79 207.41 263.88

48.72 (k ¼ π) 68.39 (k ¼π) 151.5 (k ¼ 2π) 161.0 (k ¼π)

159.5

77

symmetry group) and 8 axial terms for the fixed–free beam. For the simply supported beam, only a single axial term was used using the trigonometric approximations from the semi-analytic model. This explicitly represents the appropriate boundary conditions and reduces the problem to being purely two-dimensional. The axial response is identified by the wave number k. Six in-plane terms were used for each of the displacement variables. For both support conditions, the frequencies were computed using both the full set of elastic constants and then repeated using the elastic constants that result with the assumption of zero Poisson ratio. For the fixed–free beam, the results are shown in Table 2. As a simple comparison, the bending frequencies were also computed using the Rayleigh beam theory using the approach of Han [35]. These numerical values are the same to four significant figures with those computed by Friberg using the Vlasov beam theory [44]. The frequencies computed using the Timoshenko model of Jun [45] are in very good agreement with the values obtained using the assumption of zero Poisson ratio. There is some discrepancy with the results of Petrolo and co-workers [48], which is mildly surprising given that their approach uses similar powerseries approximations over the cross-section, albeit via forms incorporated into finite element approximations. The largest difference in frequency comes in comparison with the bending frequencies of the present model using the two constitutive theories. Those for zero Poisson ratio, which is a common assumption in many beam theories, are consistently five to ten percent lower than those computed using the true elastic stiffness tensor. The torsional frequencies are also influenced but to a much lesser extent. The results for the simply supported beam are given in Table 3, for which the main comparative study is the work of Jun [45,46]. As before, the Rayleigh theory results can be computed using the methodology of Han for the simply supported beam [35] and these results are in good agreement with the Timoshenko results of Jun and the results of de Borbon and Ambrosini [47] with somewhat less agreement with the work of Petrolo and co-workers [48]. The present model gives frequencies that are close to those of Jun but are consistently lower by up to a single percentage. Additionally, there is little change in results using the present approach when the full elastic constants are used in the calculations. Plots of the modal shapes along the length of the beam as measured at the outer diameter/thickness of the section are shown in Figs. 4–9 for the first three bending and torsional modes. The bending frequencies are those for which Rayleigh frequencies are computer and hence give modal patterns dominated by the vertical or transverse displacement component v. The torsional frequencies are dominated by the in-plane displacement components u and v but with a smaller warping motion out of the plane. 3.3. Closed sections In this section, two representative examples are studied that have been analyzed using other methods of approximation: the hollow monoclinic rectangular cantilevered beam and a homogeneous hollow elliptical beam under simple support. Both of these solids have beam-like behavior and have been studied by at least one kinematic model that restricts the displacement field. The rectangular section has relatively thick walls with anisotropic constitutive behavior. The elliptical section has isotropic properties with relatively thin walls. Hence they provide a good range of comparison with the present model that does not restrict the nature of the displacement field other than truncating the number of terms in the series solution. 3.3.1. The hollow rectangle Song and Librescu [17] considered the free vibrations of a

78

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

Fig. 3. The first four modes of the doubly unsymmetric channel. (a) Mode 1 (symmetric), (b) Mode 2 (symmetric), (c) Mode 3 (anti-symmetric), (d) Mode 4 (symmetric).

0

Table 2 Frequencies of the cantilevered thin-walled semi-circular beam.

1 2 3 4 5 6 7 8 9 10 11

Frequency [44]

Rayleigh

[45]

3DRE (ν ¼ 0)

3DRE

[48]

31.80 63.76 137.5 199.0 278.2 483.9 556.3 657.3 767.5 1075

31.80

31.8 63.79 137.7 199.3 278.4

31.77 63.66 137.1 197.5 276.6 480.7 545.9 637.9 759.7 953.2 1051.5

33.47 64.04 137.8 207.8 278.0 483.8 573.8 642.0 765.7 953.2 1103.1

31.00 67.33 357.8 193.3 593.8

199.0

556.3

558.1

1087

533.0

-0.2 Normalized displacement

Mode

-0.4

-0.6

-0.8 hollow orthotropic thin-walled beam composed of an orthotropic solid that has been rotated about a single axis and is fixed at one end. In the nomenclature of this paper, the fibers always lie in the x–y plane with the rotation angle θ being measured between the þx- and þy-axes. When θ ¼0, the fibers are perpendicular to the long axis of the beam, and when θ ¼ π/2, the fibers are then along the axis of the beam. The material properties used are those of Song and Librescu and are taken as E1 ¼206.843 GPa (30  106 psi), E2 ¼E3 ¼5.171 GPa (0.75  106), ν12 ¼ ν13 ¼ ν23 ¼0.25, G13 ¼G23 ¼2.551 GPa (0.37  106 psi), and G12 ¼3.1026 GPa (0.45  106). Here both SI and US units have been given as the latter were used in the original study [17]. The original outer dimensions

-1 0

0.2

0.4 0.6 Normalized distance

0.8

1

Fig. 4. The first symmetric mode of the semi-circle: u (short dash), v (long dash), and w (solid).

of the rectangular section were given as 0.026416 m (1.04 in) in the x-direction by 0.006096 m (0.24 in) in the y-direction with a

Table 3 Frequencies of the simply supported thin-walled semi-circular beam. Mode

1 2 3 4 5 6 7 8 9 10

Frequency λ

[45]

Rayleigh

2L 2L 2L L L 2L/3 2L/3 L/2 L 2L/5

89.27 150.4 320.3 357.1 365.8 604.1 803.5 885.0 1107 1218

89.24

356.5

800.5

[46]

[47]

[48]

3DRE (ν ¼0)

3DRE

2DFEM

149.66 317.25 364.02

89.23 149.74 317.78 356.44 364.31

86.28 181.0 278.6 784.3 353.1

89.03 149.9 315.7 352.8 364.0 599.4 781.1 874.0 1047 1195

89.06 150.0 315.9 353.1 364.6 601.3 782.9 878.2 1048 1202

89.26 149.0 314.8 353.9 361.1 596.1 784.4 871.9 1049 1195

761.8

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

0.75

79

0.25

0.5 Normalized displacement

Normalized displacement

0 0.25 0 -0.25 -0.5

-0.25

-0.5

-0.75 -0.75 -1

-1 0

0.2

0.4 0.6 Normalized distance

0.8

1

Fig. 5. The second symmetric mode of the semi-circle: u (short dash), v (long dash), and w (solid).

0

0.2

0.4 0.6 Normalized distance

0.8

1

Fig. 7. The first unsymmetric mode of the semi-circle: u (short dash), v (long dash), and w (solid).

0

0.75 0.5 Normalized displacement

-0.2

0.25 0 -0.25 -0.5

-0.4

-0.6

-0.8

-0.75 -1

-1

0

0

0.2

0.4

0.6

0.8

1

0.2

0.4 0.6 Normalized distance

0.8

1

Fig. 6. The third symmetric mode of the semi-circle: u (short dash), v (long dash), and w (solid).

Fig. 8. The second unsymmetric mode of the semi-circle: u (short dash), v (long dash), and w (solid).

thickness of 0.001016 m (0.04 in) and a length of 0.254 m (10 in) in the y-direction. The beam is fixed at y¼0. This implies that each of the displacement components is zero for all (x,z) at y¼0. The density was taken to be 1528.15 kg/m3 (14.3  10  5 lb/s2/in4). Song and Leissa incorporated a kinematic model that enforces zero deformation in the cross-sectional plane, allows for both primary and secondary warping, and assumes constant transverse shear strain. Such a model is very useful in terms of computational effort but is likely to result in an overly stiff system that overestimates the bending frequencies. Part of the objective of this work is to assess such limits. For the present model, the accuracy of the results is linked to the number of terms in the series used to represent each of the

three displacement components. As the number of terms increases, the computed frequencies get lower and more accurate until an increase in the number of terms has little influence on the results. For this geometry, an appropriate number of terms was assessed by computing the values of the lowest six frequencies for the case of θ ¼0. It was found that a combination of nine axial terms (beginning with y1 since the displacements must be zero at the fixed support) and five in-plane terms provides sufficient accuracy for all of the results that follow. For the extreme cases of θ ¼0 and θ ¼ π/2, this beam configuration has axial, torsional, and flexural modes that completely uncouple. Group theory could have been used to drastically reduce the size of the unknowns for this problem by exploiting the

80

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

1

Table 4 Axial frequencies of hollow rectangular section.

0.75

Ply angle

Mode type

Normalized displacement

First axial

Second axial

0.5 0.25

0 15 30 45 60 75 90

0 -0.25 -0.5

[17]

3DRE

1D

[17]

3DRE

11 376

11 749

34 127

35 247

71 947

72 007

11 384 11 537 12 191 13 956 18 442 35 692 71 944

215 841

216 021

34 149 34 599 36 497 41 598 53 785 100 950 216 134

Table 5 Torsional frequencies of hollow rectangular section.

-0.75

Ply angle

-1

Mode type First mode

0

0.2

0.4 0.6 Normalized distance

0.8

1

Fig. 9. The third unsymmetric mode of the semi-circle: u (short dash), v (long dash), and w (solid).

symmetry of the resulting displacements. However, this was not used to maintain the same number of degrees of freedom for random angles of fiber orientation. The axial frequencies of a bar that correspond to zero Poisson ratio approximations and a simplified displacement field are given by [49]

ωn =

1D

(2n + 1)πc 2L

(16)

Here c is the one-dimensional longitudinal wave speed given by E/ρ . E is used to denote the elastic modulus in the direction of the longitudinal wave and L is the length of the bar. For θ ¼0, E ¼E2 and for θ ¼90, E ¼E1. Similarly, the bending frequencies in the transverse (w-dominant) and lateral (u-dominant) directions can be approximated by the Rayleigh beam theory [35]. These assume axial displacements that are linear in either z or x, respectively, and also incorporate the rotary inertia of the cross-section in the kinetic energy term. However, they do not include any shear deformation. These calculations are more involved than using a single formula, but the frequencies can be easily calculated using the procedure outlined by Han and co-workers [35]. Once again, the values of modulus of elasticity and second moment of area required by the Rayleigh theory are obvious depending on the direction of bending and can be computed for the two extremes for the angle of fiber orientation. The axial frequencies are given in Table 4 for the cases of θ ¼ 0 and 90° for the one-dimensional and Song and Librescu theories, and for 15° increments for the present model for both the first and second modes. For all angles other than 0 or 90, the axial modes are not uncoupled from the rest of the motion. However, the modal displacements are the largest of all components and are generally constant over the beam cross-section making them easy to identify. Two results are notable. First, there is a surprising four percent disagreement between the result of Song and Librescu for the first axial mode when θ ¼ 0. This mode is dominated by the elastic modulus of the matrix material, and it is odd that the theory of Song and Librescu gives such a noticeable difference. The second result of note is the strongly nonlinear relationship

0 15 30 45 60 75 90

Second mode

[17]

3DRE

4239/4640

4397.6 4620.1 5266.0 6302.1 7538.7 7397.9 4686.0

[17]

3DRE

12 717/14 595

13 141 13 795 15 665 18 629 19 598 18 055 14 371

between the axial frequency and the angle of orientation of the fibers. Changing the angle from 0 to 45° results in under 25 percent increase in frequency, but an additional 45° of rotation results in an increase in frequency over a factor of 6. Hence the action of the fibers is strongly delayed until the final 15–20° of rotation align with the axial motion. The torsional frequencies are given in Table 5. The uncoupled mode for θ ¼90 is given by Song and Librescu for two conditions associated with neglecting or enforcing warping restraint within their theory. The latter case results in higher frequencies, which are in good agreement with the present model. These results also indicate that the torsional stiffness for this material peaks near 60° for the limited increments of this analysis. Of most interest are the flexural modes since they are usually the easiest to excite and possess the lowest frequency. These modes are classified according to their dominant displacement motion using the phrases transverse bending modes and lateral bending modes. These incorporate flexural motion for which the displacements in the plane of the cross-section are generally in the z (transverse) or x (lateral) directions. Once again, for fiber angles of 0 and 90° the flexural motion uncouples from the torsional and axial modes, but for other angles there is significant rotation of the cross-section coupled with the flexural motion along the lengths of the beam. The first three transverse modes are given in Table 6 and the first two lateral modes are given in Table 7. The results of Song and Librescu are compared with the present model for 15° increments, and values for the Rayleigh theory are also given for purposes of comparison for the limiting cases of 0 and 90°. As expected, the resulting frequencies from the full equations of elasticity are lower than those of the one-dimensional beam model of Song and Librescu. The agreement for the lowest transverse mode is quite good, with the beam theory results being within six percent of the elasticity results for all fiber orientations. The second transverse mode has good agreement for the lower fiber angles but overestimate the elasticity results by over twenty percent as the fibers approach the axis of the beam.

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

81

Table 6 Transverse bending frequencies of hollow rectangular section. Ply angle

Mode type First transverse

0 15 30 45 60 75 90

Second transverse

Rayleigh

[17]

3DRE

Rayleigh

[17]

3DRE

Rayleigh

3DRE

242.4

249 253 269 318 443 765 1502

242.2 245.8 262.5 306.8 415.6 730.2 1420.6

1517.2

1558 1577 1676 1985 2786 4934 8502

1497.5 1520.3 1622.0 1888.6 2525.4 4123.1 6706.7

4239.7

4106.7 4171.3 4446.6 5149.8 6730.2 10 423 14 810

1534.7

9595.9

Table 7 Lateral bending frequencies of hollow rectangular section. Ply angle

Mode type First lateral

0 15 30 45 60 75 90

Second lateral

Rayleigh

[17]

3DRE

Rayleigh

[17]

3DRE

848.4

871 881 949 1236 2182 3916 4583

844.5 859.1 917.8 1069.7 1438.3 2484.6 4476.2

5234.5

4331 4383 4657 5490 7554 12 333 18 062

5059.7 5170.3 5560.1 6462.7 8473.3 12 944 17 269

5365.6

Third transverse

33 106

The lateral frequencies are in good agreement for fiber orientations near 0 or 90°, but for other orientations are as much as 50 percent higher than the elasticity results. The second lateral modes are more problematic. Both the Rayleigh theories and elasticity theory give a frequency of over 5000 rad/s for this mode when the fibers are oriented at 0°. The results of Song and Librescu have a frequency of 4331 rad/s. This pattern continues for all fiber orientations in that the results of Song and Librescu are roughly

26 814

10–15 percent under the elasticity results. This surprising result may be explained if the second lateral mode identified by Song and Librescu was actually the third transverse bending mode, in which case the comparison is far more in line with the elasticity results given in Tables 6 and 7. Typical modal shapes for two representative fiber orientations are shown in Figs. 10–13. The lowest two torsional modes are shown in Fig. 10, and demonstrate the coupled nature of the deformation as the fiber angle goes from 0° to 60°. In the former case, the torsional motion is dominant and virtually exclusive. As the fiber angle changes, the torsional motion still dominates the nature of the cross-sectional deformation but there is far more inclusion of deformation that is not strictly torsional that contributes to the stored energy and explains in part the increase in frequency for orientation angles between 0 and 90°. Similar behavior is demonstrated in Figs. 11 and 12, which show the lowest two transverse bending modes. For θ ¼0, the motion is purely flexural. For other fiber angles, there is combined visible twist of the crosssection along with the dominant bending action. In the case of bending, the fibers are being reoriented to change what amounts to the longitudinal modulus of the beam, which also increases the frequency. Only in the case of the lowest two lateral modes shown in Fig. 13 do the modes stay visually uncoupled. This reflects the

Fig. 10. The first (top) and second (bottom) torsional modes for two fiber orientations. (a) θ ¼ 0, (b) θ ¼ 60°.

82

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

Fig. 11. The first transverse bending mode from the side (top) and end (bottom) for two fiber orientations. (a) θ ¼ 0, (b) θ ¼ 60°.

fact that the fibers are in the plane of the primary cross-section deformation rather than perpendicular to the motion as is the case of transverse bending. The original material properties of Song and Librescu were somewhat unusual for a fiber-reinforced material in the Poisson ratios were identical and the shear modulii in the plane of the fibers were unequal. To present additional results for another material, the analysis with the Song and Librescu geometry was repeated for the graphite–polymer material properties tabulated by Hyer [50]: E1 ¼155.0 (all GPa), E2 ¼E3 ¼12.10, ν12 ¼ ν13 ¼0.248, ν23 ¼ 0.458, G13 ¼G12 ¼4.40, G23 ¼ 3.20. A density of ρ ¼1600 kg/m3 was used, which is typical for this class of material. The results of this analysis are given in Table 8, with the values given for the first six frequencies for fiber orientation angles of π/12. These are mainly given for purposes of comparison and to serve as a benchmark for other computational models. One crucial difference between the results using the material properties of Hyer and the general behavior of Song and Librescu is the nature of the change in frequency of the primary torsion mode as the fiber orientation angle changes. For the material properties of Song and Librescu, this frequency peaked at θ ¼ 45° with an increase in frequency of about 25 percent. For the properties of Hyer, this mode peaks at about 66° with a frequency of 8780 rad/s. This is a 73 percent increase over the value at θ ¼0. When the fibers align with the global x or z coordinates (θ ¼0 or θ ¼ 90, respectively), the frequency depends almost exclusively on the shear modulii in the 23 plane (a lower value) and then the 13 plane (a higher value). 3.3.2. The hollow ellipse The vibration of hollow shells of elliptic cross section, often referred to as elliptic-cylindrical (EC) shells when the wall thickness is thin, has a number of practical applications even though the number of comparative studies is relatively small. For shells with constant thickness, Sewall and co-workers [51] provided both experimental and analytical results and Yamada and co-workers [52] gave numerical results using classical shell theory. More

recently Hayek and Boisvert [53] studied similar EC shells with a higher-order shell theory. The shell geometry and constitutive properties considered in this section matches that of an identical shell studied by both Yamada and co-workers [52] and by Hayek and Boisvert [53]. The ratio of the semi-major axis length a to that of the semi-minor axis length of b is given as a/b = 2.5, and the ratio L/a = 2.198, where L is the length of the beam. The ratio of the shell thickness h to a is given by h/a = 0.007325. Both ends of the cylinder are under simple support. Two cases are considered: the isotropic shell and the orthotropic shell. A general spatial solution to this class of motion can be expressed as ∞

u(x, y , z ) =

∑

Um(x, z )sin

mπy L

Vm(x, z )cos

mπy L

m=1 ∞

v(x, y , z ) =

∑ m=0 ∞

w (x, y , z ) =

∑

Wm(x, z )sin

m=1

mπy L

These displacements satisfy the boundary conditions at both ends of the shell and, unlike the cantilever section, split the axial and cross-sectional dependencies of the problem. This allows frequencies to be determined for specified values of m. For this specific problem, the dimensionless frequencies Ω are sought for the isotropic case when the Poisson ratio ν ¼ 0.3 and m ¼ 1. They are computed as dimensionless frequencies, and are given by [53]

Ω2 =

ω2a2ρ(1 − ν 2) E

(17)

where E is the modulus of elasticity. As was the case for the hollow rectangle, the numerical integration over the thin walled section can be accomplished in closed form by evaluating the integrals of the coefficient matrices for the outer semi-axis lengths and then subtracting similar integrals over the inner semi-axis lengths.

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

83

Fig. 12. The second transverse bending mode from the side (top) and end (bottom) for two fiber orientations. (a) θ ¼ 0, (b) θ ¼ 60°.

Since this section is isotropic, the vibrational modes can be separated into various modal groups that isolate the nature of the motion. In the studies of others, the displacement functions over the beam cross-section were explicitly separated in terms of trigonometric functions that describe the circumferential variation around the circumferential coordinate θ. This effectively makes the problem one-dimensional in the circumferential direction of the cross-section mid-plane. Both Yamada [52] and Hayek and Boisvert [53] reported their results as a function of either symmetric or antisymmetric modes along with the circumferential wave number n in the above expansion. In the present analysis, the frequencies are grouped according to the symmetries of the three displacement components (u,v,w) in the cross-section (x,z) directions. A summary of these groups is given in Table 9, which enforces antisymmetry about an axis by the symbol O, indicating only odd polynomials for that displacement component were given in that coordinate direction. The numerical results are presented in terms of these four groups. As was the case for the rectangular section, the number of terms used to represent the displacement components was

adjusted until an increase in number of terms did not yield a significant drop in frequency. Since there is now a single axial term containing either the sine or cosine axial dependence, only the cross-section variables need adjustment. Approximations up to and including x9z7 for u, x9z8 for w, and x8z9 for v were used in the results that follow. The dimensionless frequencies for all four groups are tabulated in Table 10 and compared with the results of Yamada [52] and Hayek and Boisvert [53] along with the two-dimensional finite element model developed as part of the present study. These latter results were obtained by meshing one-quarter of the elliptical cross-section with 4 divisions in the thickness direction and 541 elements along one quadrant of the ellipse. These numbers were determined after both halving and doubling the number of divisions with little change in frequency. Hence there were a total of 2160 4-node elements and 2705 nodes used to represent the quadrant. Four different analyses were ran, each of which had boundary conditions at the 0 and 90° borders of the quadrant that matched those of the grouped polynomials. The results are all in very good agreement with the results from

84

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

Fig. 13. The lowest lateral bending modes from the side (top) and end (bottom). These modal plots are visually identical for both θ¼ 0 and θ¼ 60. (a) First mode, (b) second mode.

Table 8 The first six modes of hollow rectangle with modified elastic constants. Angle

1

2

3

0 (Rayleigh)

362.8

1268

2268

0 15 30 45 60 75 90

363.1 356.9 354.8 385.8 489.0 779.4 1252

1264 1244 1236 1337 1686 2658 4094

2223 2191 2187 2379 2992 4517 5460

90 (Rayleigh)

1298

4539

8118

4

5086 5459 5980 6519 8093 7980 6649

Table 10 Frequencies of hollow isotropic beam with elliptical cross-section.

5

6

Ply angle

Method

6338

7825

Group

3DRE

[53]

[52]

2d FEM

6013 5951 6471 7861 8982 11720 15702

7431 7388 7474 8110 10001 14202 16697

1 1 1

0.0706 0.1285 0.1631

0.0722 0.1327 0.1644

0.0721 0.1327 0.1644

0.0719 0.1319 0.1641

2 2 2

0.0762 0.1098 0.1810

0.0778 0.1111 0.1877

0.0782 0.1111 0.1874

0.0776 0.1118 0.1872

22 685

28 008

3 3 3

0.0763 0.1087 0.1793

0.0783 0.1119 0.1872

0.0778 0.1110 0.1879

0.0777 0.1107 0.1856

4 4 4

0.0702 0.1291 0.2026

0.0722 0.1340 0.2012

0.0721 0.1339 0.2005

0.0716 0.1332 0.2023

Table 9 Modal symmetry groups for beam with elliptical cross-section. Group

Displacement

x

u v w

O

1

2

u v w

3

u v w

4

u v w

z

O

O O

O

O

O O O

O O

O

theory for the third modes being closer to 4 percent difference. This good agreement is to be expected since the shell in this case is quite thin and hence its kinematics should be captured by such a theory. Hence when the beam is isotropic, the simplified shell theory overpredicts the elasticity frequencies by approximately three percent. For orthotropic materials, the analysis was repeated using the material properties of Hyer [50] listed for the rectangular section. The fibers are oriented along the beam axis (θ ¼90). In this case the same symmetry groups hold as for the isotropic solid. The dimensionless frequencies in this case are given by

Ω2 = the elasticity model giving, with only one exception, the lower frequency. The average differences between the lowest 3 frequencies from each symmetry group is 2.5, 2.8, and 2.1 percent between the elasticity results and the higher-order shell theory of Hayek and Boisvert [53] with some of the frequencies from shell

ω2a2ρ E1

(18)

These values are shown in Table 11 and are again compared with the results from the semi-analytic finite element model. These frequencies are significantly below those of the isotropic section. This is partially influenced by the ratio of the longitudinal modulus

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

For closed thin-walled sections:

Table 11 Frequencies of hollow orthotropic beam with elliptical cross-section. Group

85

Method 3DRE

2D FEM

1 1 1

0.0288 0.0419 0.0619

0.0292 0.0428 0.0631

2 2 2

0.0286 0.0446 0.0595

0.0290 0.0455 0.0604

3 3 3

0.0285 0.0438 0.0564

0.0289 0.0445 0.0574

4 4 4

0.0289 0.0419 0.0661

0.0293 0.0431 0.0665

to that of the shear modulus. For the isotropic solid, this ratio is 2.6. For the orthotropic material, the ratio is over 35. Hence the level of shear deformation is expected to be significantly higher than that of the isotropic solid, and any simplified model would need to take this behavior into account.

4. Conclusions An elasticity-based approximation scheme was used to develop estimates for the frequencies of a wide variety of thin-walled beams. In nearly every case, this method provided lower estimates of all frequencies studied. The primary conclusions of the application of this approach as applied to open and closed sections are given below. For open thin-walled sections:

1. For a hollow anisotropic rectangular section whose thickness to cross-section length direction was 1:6 for one wall and 1:21 for another, good agreement was found for axial and torsional modes when fibers had 0 or 90° orientation angles between elasticity theory and a refined beam theory. 2. The lowest (transverse) bending frequencies of the rectangular section were also in good agreement for all fiber orientations, with a maximum difference of about 7 percent. The second transverse frequencies had more significant differences that increased with fiber orientation, starting at 4 percent (θ ¼0) and increasing to over 25 percent (θ ¼ 90). 3. The lateral bending frequencies of the rectangular section were in good (under 4 percent difference) agreement for the cases where the fibers were along or perpendicular to the beam axis, but for other angles of orientation the present method gave frequencies that were significantly (over 50 percent) lower than those of the one-dimensional beam model. It is possible that the second lateral mode of Song and Librescu [17] is actually the third transverse mode. 4. For the isotropic elliptical shell whose thickness to cross-section dimension ratio was 1:136, excellent agreement was obtained with thin-shell theory results with most differences in frequency being less than 3 percent. 5. The dimensionless frequencies of an orthotropic elliptical beam are roughly a third of those for the isotropic beam, indicating a dramatic loss in general stiffness even when the wall thickness is small. In addition to these conclusions, representative elasticity results were presented for anisotropic (for the rectangle) and orthotropic (for the ellipse) thin-walled beams that can be used for comparison with other simplified models.

Acknowledgment 1. The frequency results of the present model for the doubly unsymmetric section of Yaman [39] are in very good agreement with two and three-dimensional finite element results and the lowest modes of other researchers but are significantly lower for the remaining modes. 2. For the cantilevered thin-walled semi-circle with an assumed Poisson ratio of zero, the bending frequency results computed using the elasticity-based model are similar to but slightly lower than results from Vlasov and Timoshenko beam theories, the former of which are exactly the same for the beam studied as the results from the Rayleigh theory. A similar pattern is true for torsional modes. 3. For the cantilevered thin-walled semi-circle with a non-zero Poisson ratio, the bending frequency results from the elasticity model are slightly above (roughly five percent) those computed by Vlasov and Timoshenko beam theories. There is little difference between torsional frequencies for the present model in terms of influence of Poisson ratio. 4. For the simply supported thin-walled semi-circle, the results of the elasticity model are slightly below those of Vlasov and Timoshenko model results regardless of the value of Poisson ratio used. The assumption of zero Poisson ratio in the elasticity model developed here has a very small (well under one percent) influence on the simply supported frequencies. 5. The agreement between the present elasticity model used here and a similar power series approach presented by Petrolo and co-workers [48] is mixed, and may have been influenced by either truncated approximations or the influence of local axial polynomials rather than the global polynomials used in the present study.

This work was sponsored by the Mountains-Plains Consortium. The support is gratefully acknowledged. This paper is dedicated to the courageous example of Professor Liviu Librescu.

Appendix The elements of the matrix equations for both the Ritz model and the three-dimensional finite element model are given as

[K11]ij =

⎡

u

∫V ⎢⎢⎣C11∂∂Ψxi

∂Ψ ju ∂x

u

+ C16

u

∂Ψiu ∂Ψ j ∂Ψ u ∂Ψ j + C55 i ∂x ∂y ∂z ∂z

u u ∂Ψi ∂Ψ j ∂Ψ u ∂Ψ j ⎤ ⎥ dV + C66 i ∂y ∂x ∂y ∂y ⎥⎦ u

+ C16

[K12]ij =

⎡

u

∫V ⎢⎢⎣C12 ∂∂Ψxi

∂Ψ jv ∂y

v

+ C16

(19) v

∂Ψiu ∂Ψ j ∂Ψ u ∂Ψ j + C45 i ∂x ∂x ∂z ∂z

v v ∂Ψi ∂Ψ j ∂Ψ u ∂Ψ j ⎤ ⎥ dV + C66 i ∂y ∂y ∂y ∂x ⎥⎦ u

+ C26

[K13]ij =

⎡

u

∫V ⎢⎢⎣C13 ∂∂Ψxi + C36

∂Ψjw ∂z

w

+ C45

(20) w

∂Ψiu ∂Ψj ∂Ψ u ∂Ψj + C55 i ∂z ∂y ∂z ∂x

w⎤

∂Ψiu ∂Ψj ∂y ∂z

⎥ dV ⎥⎦

(21)

86

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

[K22]ij =

⎡

v

∫V ⎢⎢⎣C22 ∂∂Ψyi

∂Ψ jv ∂y

v

+ C26

v

∂Ψiv ∂Ψ j ∂Ψ v ∂Ψ j + C44 i ∂y ∂x ∂z ∂z

v v ∂Ψ ∂Ψ j ∂Ψ v ∂Ψ j ⎤ ⎥ dV + C26 i + C66 i ∂x ∂y ∂x ∂x ⎥⎦

[3]

v

(22)

[4] [5]

w w ⎡ ∂Ψ v ∂Ψ w ∂Ψ v ∂Ψj ∂Ψ v ∂Ψj j ⎢C23 i [K ]ij = + C44 i + C45 i V ⎢ ∂y ∂z ∂z ∂y ∂z ∂x ⎣ 23

∫

+ C36

w ∂Ψiv ∂Ψj ⎤ ⎥ dV ∂x ∂z ⎥⎦

[6] [7]

(23)

[8] [9]

w w ⎡ ∂Ψ w ∂Ψ w ∂Ψ w ∂Ψj ∂Ψ w ∂Ψj j ⎢C33 i [K 33]ij = + C44 i + C45 i V ⎢ ∂z ∂z ∂y ∂y ∂y ∂x ⎣

[10]

∫

w w ∂Ψ w ∂Ψj ∂Ψ w ∂Ψj ⎤ ⎥ dV + C45 i + C55 i ∂x ∂y ∂x ∂x ⎥⎦

[11] [12]

(24) [13]

Mij11

∫V

=

ρϕi uϕju

dV

(25)

[14] [15]

Mij22 =

∫V ρϕivϕjv dV

Mij33 =

∫V ρϕi wϕjw dV

(26)

[16] [17]

(27)

[18] [19]

The elements of the coefficient matrices for the two-dimensional semi-analytic finite element model are given as

Kij11 =

⎛

∫V ⎜⎜C11∂∂Φxi

u

⎝

+ C55

⎛

Kij13 =

∫V ⎜⎜⎝C13 ∂xi

⎛

[22]

cos2ky + C66ΦiuΦjusin2ky

[23]

(28)

⎞ ∂Φjv ∂Φiu v ksin2ky⎟⎟ dV = K21 Φj kcos2ky − C66Φiu ji ∂x ∂x ⎠

∫V ⎜⎜⎝C12

w

(29)

(30)

v ⎞ ∂Φiv ∂Φj sin2ky⎟⎟ dV . ∂x ∂x ⎠

⎛

∂Φ v

Φjwkcos2ky − C44Φiv

[30] [31]

(31)

⎞ ksin2ky⎟⎟dV = K 32 ji ∂z ⎠

(32)

[34] [35]

∫

w ⎞ ∂Φi w ∂Φj cos2ky+⎟⎟ dV ∂x ∂x ⎠

[32]

[33]

∂Φjw

⎛ ∂Φ w ∂Φjw ⎜⎜C33 i cos2ky + C44Φi wΦjwk 2sin2ky = V ⎝ ∂z ∂z + C55

[28] [29]

∫

∫V ⎜⎜⎝C23 ∂zi

[26]

[27]

⎟⎟cos2ky dV = K 31 ji ⎠

v ⎛ ∂Φ v ∂Φj ⎜⎜C22ΦivΦjvk 2cos2ky + C44 i sin2ky = V ⎝ ∂z ∂z

Kij23 =

[24] [25]

w⎞

∂Φ u ∂Φj ∂Φ u ∂Φj + C55 i ∂z ∂z ∂x

+ C66

Kij33

∂x

u ⎞ ∂Φiu ∂Φj cos2ky⎟⎟ dV ∂z ∂z ⎠

Kij12 =

Kij22

∂Φju

[20] [21]

[36] [37]

(33)

[38] [39]

References [1] H. Wagner, W. Pretschner, Verdrehung und Knickung von offenen Profilen, Luftfahrtforschung 11 (1934) 174–180. [2] F. Bleich, H. Bleich, Bending, torsion, and buckling of bars composed of thin

[40]

[41]

walls, Prelim. Pub. 2nd Congress of International Association for Bridge and Structural Engineers, English edition, 1936, p. 871. R. Kappus, Twisting failure of centrally loaded open-section columns in the elastic range. Luftfahrtforschung 14 (1937) 444–457. (Translated by J. Vanier, National Advisory Committee for Aeronautics, No. 851 (1938)). A. Gjelsvik, The Theory of Thin-walled Bars, John Wiley and Sons, New York, 1981. V. Vlasov, Thin-walled Elastic Beams, 2nd ed., Israel Program for Scientific Translations, Jerusalem, 1963. N.J. Goodier, The Buckling of Compressed Bars by Torsion and Flexure, Cornell University Engineering Experiment Station, Bulletin 27, 1941. S.P. Timoshenko, Theory of bending, torsion and buckling of thin-walled members of open cross section, J. Frankl. Inst. 239 (1945) 201–219. J.M. Gere, Y. Lin, Coupled vibrations of thin-walled beams of open cross-section, J. Appl. Mech. ASME 25 (1958) 373–378. J.R. Banerjee, S. Guo, W.P. Howson, Exact dynamic stiffness matrix of a bendingtorsion coupled beam including warping, Comput. Struct. 59 (1996) 613–621. C. Mei, Coupled vibrations of thin-walled beams of open section using finite element method, Int. J. Mech. Sci. 12 (2001) 883–891. L.P. Kollar, Flexural–torsional vibration of open section composite beams with shear deformation, Int. J. Solids Struct. 38 (2001) 42–43. V.H. Cortinez, M.T. Piovan, Vibration and buckling of composite thin-walled beams with shear deformability, J. Sound Vib. 258 (2002) 701–723. M.T. Piovan, V.H. Cortinez, Mechanics of shear-deformable thin-walled beams made of composite materials, Thin-Walled Struct. 45 (2007) 37–62. S.P. Machado, V.H. Cortinez, Free vibration of thin-walled composite beams with static initial stresses and deformations, Eng. Struct. 29 (2007) 372–382. S.U. Benscoter, A theory of torsion bending for multicell beams, ASME J. Appl. Mech. 21 (1954) 25–34. D.H. Hodges, A.R. Atilgan, M.V. Fulton, L.W. Rehfield, Free-vibration analysis of composite beams, J. Am. Helicopter Soc. 36 (1991) 36–47. O. Song, L. Librescu, Free vibration of anisotropic composite thin-walled beams of closed cross-section contour, J. Sound Vib. 167 (1993) 129–147. N.W. Murray, Introduction to the Theory of Thin-walled Structures, Clarendon Press, Oxford, 1984. L. Librescu, O. Song, Thin-walled Composite Beams: Theory and Application, Springer, The Netherlands, 2006. D.H. Hodges, Nonlinear Composite Beam Theory, AIAA, United States, 2006. M. Hajianmaleki, M.S. Qatu, Vibrations of straight and curved composite beams: a review, Compos. Struct. 100 (2013) 218–232. W.M. Visscher, A. Migliori, T.M. Bell, R.A. Reinert, On normal modes of free vibration of inhomogeneous and anisotropic elastic objects, J. Acoust. Soc. Am. 90 (1991) 2154–2162. J.N. Reddy, Energy Principles and Variational Methods in Mechanics, Wiley, New York, 2002. J.N. Reddy, Mechanics of Laminated Composite Plates and Shells: Theory and Analysis, CRC Press, Boca Raton, 2004. H.H. Demarest, Cube resonance method to determine elastic constants of solids, J. Acoust. Soc. Am. 49 (1971) 768–775. I. Ohno, Free vibration of a rectangular parallelepiped crystal and its application to determination of elastic constants of orthorhombic crystals, J. Phys. Earth 24 (1976) 355–379. P.R. Heyliger, Axisymmetric free vibrations of finite anisotropic cylinders, J. Sound Vib. 148 (1991) 507–520. P.R. Heyliger, A. Jilani, The free vibrations of inhomogeneous cylinders and spheres, Int. J. Solids Struct. 29 (1992) 2689–2708. A. Migliori, J.L. Sarrao, Resonant Ultrasound Spectroscopy, Wiley, New York, 1997. P. Heyliger, P. Ugander, H. Ledbetter, Anisotropic elastic constants: measurement by impact resonance, ASCE J. Mater. Civil Eng. 13 (2001) 356–363. P.R. Heyliger, Ritz finite elements for curvilinear particles, Commun. Numer. Methods Eng. 22 (2005) 335–345. E. Carrera, M. Petrolo, P. Nali, Unified formulation applied to free vibrations finite element analysis of beams with arbitrary cross-section, Shock Vib. 18 (2010) 485–502. E. Carrera, G. Giunta, M. Petrolo, Beam Structures: Classical and Advanced Theories, John Wiley and Sons, New York, 2011. E. Mochizuki, Application of group theory to free oscillations of an anisotropic rectangular parallelepiped, J. Phys. Earth 35 (1987) 159–170. S.M. Han, H. Benaroya, T. Wei, Dynamics of transversely vibrating beams using four engineering beam theories, J. Sound Vib. 225 (1999) 935–988. S.P. Timoshenko, On the transverse vibrations of bars of uniform cross-section, Philos. Mag. 125 (1922) 125–131. S.P. Timoshenko, On the correction for shear of the differential equation for transverse vibrations of bars of uniform cross-section, Philos. Mag. 41 (1921) 744–746. G.R. Cowper, The shear coefficient in Timoshenko's beam theory, J. Appl. Mech. 33 (1966) 335–340. Y. Yaman, Vibrations of open-section channels: a coupled flexural and torsional wave analysis, J. Sound Vib. 204 (1997) 131–158. Y. Yaman, personal communication (2014). The original dimension of the thickness was 0.05 inches and the centerline leg lengths of the channel were 0.5, 1.5, and 1 inches. These values are extremely close to the full dimensions of the cross-section used for purposes of comparison in this study. D. Ambrosini, On free vibration of nonsymmetrical thin-walled beams, ThinWalled Struct. 47 (2009) 629–636.

P.R. Heyliger / Thin-Walled Structures 95 (2015) 73–87

[42] A. Arpaci, E. Bozdag, On free vibration analysis of thin-walled beams with nonsymmetrical open cross-sections, Comput. Struct. 80 (2002) 691–695. [43] M. Tanaka, A. Bercin, Free vibration solution for uniform beams of nonsymmetrical cross-section using Mathematica, Comput. Struct. 71 (1999) 1–8. [44] P.O. Friberg, Beam element matrices derived from Vlasov's theory of open thin-walled elastic beams, Int. J. Numer. Methods Eng. 21 (1985) 1205–1228. [45] L. Jun, L. Wanyou, S. Rongying, H. Hongxing, Coupled bending and torsional vibration of nonsymmetrical axially loaded thin-walled Bernoulli–Euler beams, Mech. Res. Commun. 31 (2004) 697–711. [46] L. Jun, S. Rongying, X. Jin, Coupled bending and torsional vibration of nonsymmetrical axially loaded thin-walled Timoshenko beams, Int. J. Mech. Sci. 46 (2004) 299–320. [47] F. de Borbon, D. Ambrosini, On free vibration analysis of thin-walled beams axially loaded, Thin-Walled Struct. 48 (2010) 915–920.

87

[48] M. Petrolo, E. Zappino, E. Carrera, Refined free vibration analysis of one-dimensional structures with compact and bridge-like cross-sections, ThinWalled Struct. 56 (2012) 49–61. [49] S.S. Rao, Vibration of Continuous Systems, Wiley, New Jersey, 2007. [50] M.W. Hyer, Stress Analysis of Fiber-Reinforced Composite Materials, McGrawHill, Boston, 1998. [51] J.L. Sewall, W.M. Thompson Jr., C.G. Pusey, An experimental and analytical vibration study of elliptical cylindrical shells. NASA TN D-6089, 1971. [52] G. Yamada, T. Irie, S. Notoya, Natural frequencies of elliptical cylindrical shells, J. Sound Vib. 101 (1985) 133–139. [53] S.I. Hayek, J.E. Boisvert, Vibration of elliptic cylindrical shells: higher order shell theory, J. Acoust. Soc. Am. 128 (2010) 1063–1072.