Research Article  Open Access
HigherOrder Hierarchical Models for the Free Vibration Analysis of ThinWalled Beams
Abstract
This paper addresses a free vibration analysis of thinwalled isotropic beams via higherorder refined theories. The unknown kinematic variables are approximated along the beam cross section as a order polynomial expansion, where is a free parameter of the formulation. The governing equations are derived via the dynamic version of the Principle of Virtual Displacements and are written in a unified form in terms of a “fundamental nucleus.” This latter does not depend upon order of expansion of the theory over the cross section. Analyses are carried out through a closed form, Naviertype solution. Simply supported, slender, and short beams are investigated. Besides “classical” modes (such as bending and torsion), several higher modes are investigated. Results are assessed toward threedimensional finite element solutions. The numerical investigation shows that the proposed Unified Formulation yields accurate results as long as the appropriate approximation order is considered. The accuracy of the solution depends upon the geometrical parameters of the beam.
1. Introduction
Many typical aeronautical and space structures involve lightweight, thinwalled, beamlike structures that operate in complex environments. The free vibration characteristics are of fundamental importance in the design of such structures. In particular, the dynamic behaviour of thinwalled beams is richer than that of solid prismatic beams since, besides classical modes (such as bending, torsion, and axial deformation), local higher modes are present. Furthermore, the number of higher modes occurring before the axial one increases as the lengthtoside ratio decreases. For these reasons, the accurate modelling and analysis of thinwalled beams (as proposed in this paper) represent an interesting and uptodate research topic.
A brief overview of recent works about the free vibration of thinwalled beams follows. Matsunaga [1] analysed the natural frequencies and buckling loads of simply supported beams subjected to initial axial forces. Thin rectangular cross sections were investigated. A bidimensional displacement field was assumed. Chen et al. [2] combined the state space method with the differential quadrature method to obtain a semianalytical method for the free vibration analysis of straight isotropic and orthotropic beams with rectangular cross sections. A discussion about properties of the natural frequencies and modes for a Timoshenko beam was presented by van Rensburg and van der Merwe [3]. Attarnejad et al. [4] introduced the basic displacement functions that are calculated solving the governing differential equations of transverse motion of Timoshenko’s beams by means of the power series method. These functions were applied to the free vibration analysis of nonprismatic beams. Gunda et al. [5] analysed the large amplitude free vibration of Timoshenko beams using a finite element formulation. Transverse shear and rotatory inertia were both accounted for. Different boundary conditions were investigated. Benamar et al. [6] presented a general model for large vibration of thin straight beams. Hamilton’s principle was used to obtain a set of nonlinear algebraic equations. Simply supported and clampedclamped boundary conditions were investigated. Tanaka and Bercin [7] studied the natural frequencies of uniform thinwalled beams without cross sectional symmetry. A study on the free vibration of axially loaded slender thinwalled beams was presented by Jun et al. [8]. The effects of the warping stiffness and the axial force were included within EulerBernoulli’s beam theory. Chen and Hsiao [9] investigated these axial and torsional vibration modes of beams with a Zshaped cross section. The governing equations were derived by the principle of virtual work. The same authors presented in [10] a finite element formulation for the coupled free vibration analysis of thinwalled beams with generic open cross sections. Duan [11] studied the nonlinear free vibrations of thinwalled curved beams with unsymmetric open cross sections via the finite element method. Vörös [12] accounted for the coupling between different vibration modes due to the eccentricity of the cross section as well as to steadystate lateral loads and internal stress resultants. The governing differential equations and boundary conditions were derived using a linearised theory of large rotations and small strains and the principle of virtual work. A finite element model with seven degrees of freedom per node was formulated. Ambrosini [13, 14] presented a numerical and experimental study on the natural frequencies of beams with unsymmetric thinwalled open cross sections. Vlasov’s theory was modified to include the effects of shear deformation, rotatory inertia, and variable cross section properties. An extension of the previous theory was proposed by de Borbón and Ambrosini [15] to investigate the natural frequencies of thinwalled beams under axial loads. Experimental tests were also carried out in order to verify the proposed theory. Arpaci et al. [16, 17] investigated the free vibrations of nonsymmetric thinwalled beams with open cross sections by solving exactly the bending and torsional dynamic equilibrium equations. The coupled bendingtorsional behaviour of beams was studied by Banerjee [18] via the dynamic stiffness matrix method. S.B. Kim and M.Y. Kim [19] carried out a free vibration as well as a stability analysis of thinwalled tapered beams and frames by means of the finite element method. ZhouLian et al. [20] analysed the free vibrations in orthotropic membranes accounting for large deflections. A similar investigation was carried out by by Liu et al. [21] using the LP perturbation method.
A free vibration analysis of thinwalled isotropic beams is addressed within this paper. The relevance of the proposed analysis is due to the fact that the thinness of the cross section elements results (besides classical bending, torsion, and axial deformation) in local higher modes and corresponding frequencies. The dynamics of thinwalled beams is richer than that of solid prismatic beams and the number of higher modes before the axial one increases as the lengthtoside ratio decreases. Therefore, accurate higherorder theories are required. Furthermore, the considered analysis represents a severe test for the proposed models in order to outline both their merits and limitations. In this work, several hierarchical models as well as the classical theories are derived via a Unified Formulation (UF) that has been previously formulated for plates and shells (see Carrera [22] and Carrera and Giunta [23, 24]) and lately extended to beams with solid and thinwalled cross sections; see Carrera and Giunta [25], Carrera et al. [26], and Giunta et al. [27]. Within this UF, the a priori approximation of the displacement field is written in a compact form. The governing equations are derived through the Principle of Virtual Displacements in terms of a “fundamental nucleus.” This nucleus is an invariant of the formulation in the sense that it does not depend upon the theory order of expansion. As a result, a very broad set of beam models that account for transverse shear deformation and cross section in and outofplane warping (although a warping function is not explicitly assumed) can be obtained. Governing differential equations are solved via a Naviertype closed form solution. Open and closed thinwalled simply supported beams are investigated. Slender as well as short simply supported beams are investigated. Open and closed thinwalled cross sections are considered. The modes up to the axial one (which is not included for the sake of brevity) are all considered. The axial mode is disregarded since it is accurately predicted by classical EulerBernoulli’s and Timoshenko’s models. Results are compared with threedimensional finite element models showing that the dynamic behaviour of thinwalled beam structures can be accurately predicted.
2. Preliminaries
A beam is a structure whose axial extension () is predominant if compared to any other dimension orthogonal to it. The cross section () is identified by intersecting the beam with planes that are orthogonal to its axis. A Cartesian reference system is adopted:  and axes are two orthogonal directions laying on . The coordinate is coincident to the axis of the beam. Cross sections that are obtained by the union of rectangular subdomains withare considered; see Figure 1. Points are the coordinates of the corner points of a subdomain. Through the paper, superscript “” represents a cross section subdomain index. The cross section is considered to be constant along . The displacement field is in which , , and are the displacement components along , , and axis, respectively, and is the time variable. Superscript “” represents the transposition operator. Stress, , and strain, , vectors are grouped into vectors , that lay on the cross section: and , laying on planes orthogonal to : The straindisplacement geometrical relations are Subscripts “,” “,” and “,” when preceded by comma, represent derivation versus the corresponding spatial coordinate. A compact vectorial notation can be adopted for (6): where , , and are the following differential matrix operators: being the unit matrix.
(a)
(b)
Under the hypothesis of linear elastic materials, the generalised Hooke law, holds. The nonzero terms of the material stiffness matrix are being the Young modulus, the Poisson ratio, and the shear modulus. According to stresses and strain arrangements in (4) and (5), (9) is rewritten as follows: where matrices , , , and are
3. Hierarchical Beam Theories
The variation of the displacement field over the cross section can be postulated a priori. Several displacementbased theories can be formulated on the basis of the following generic kinematic field: are unknown functions of the spatial coordinate and time . stands for the number of unknowns that depends on the approximation order (this latter is a free parameter of the formulation). The compact expression in (13) is based upon Einstein’s notation: an index that is repeated twice implicitly stands for summation. Thanks to this notation, problem’s governing differential equations and boundary conditions can be derived in terms of a single “fundamental nucleus.” The theoretical formulation is valid for the generic approximation order and approximating functions . In this paper, the functions are the Mac Laurin polynomials. This choice is inspired by the classical beam models. and as functions of can be obtained via Pascal’s triangle as shown in Table 1. The actual governing differential equations and boundary conditions due to a fixed approximation order and polynomials type are obtained straightforwardly via summation of the nucleus corresponding to each term of the expansion. According to the previous choice for the polynomial functions, a generic order displacement field is As far as the firstorder approximation is concerned, the kinematic field is Classical models, such as Timoshenko’s beam theory (TBT), and EulerBernoulli beam theory (EBT), are straightforwardly derived from the firstorder approximation model. Within this paper no shear correction coefficient is considered for TBT, since it depends upon several parameters, such as the geometry of the cross section (see Cowper [28] and Murty [29]). Higherorder models yield a more detailed description of the shear mechanics (no shear correction coefficient is required), of the in and outofsection deformations, of the coupling of the spatial directions due to Poisson’s effect, and of the torsional mechanics than classical models do. EBT theory neglects them all since it was formulated to describe the bending mechanics. TBT model accounts for constant shear stress and strain components. It should be noted that, for TBT and EBT, the reduced Hooke law for the axial stress/strain relation should be used in accordance with the relaxation procedure classically used to contrast Poisson’s locking; see Carrera and Brischetto [30].

4. Governing Equations
The dynamic version of the Principle of Virtual Displacements (also known as d’Alembert principle) reads where stands for a virtual variation, represents the strain energy, and is the inertial work.
4.1. Virtual Variation of the Strain Energy
According to the grouping of the stress and strain components in (4) and (5), the virtual variation of the strain energy is considered as the sum of two contributions: By substitution of the geometrical relations, (7), the material constitutive equations, (11), and the unified hierarchical approximation of the displacements, (13), and after integration by parts, (20) reads In a compact vectorial form,The components of the differential stiffness matrix are The generic term is a cross section moment: that is, a weighted sum (in the continuum) of each elemental cross section area. The weight accounts for material and geometrical distribution. The integral in (24) can be easily computed for the considered domains; see (1) and (2). The following general form for the product can be assumed: where , , , and are constant depending upon indexes and as in Table 1 and whether differentiation with respect to and should be performed or not. The closed form of the generic higherorder moment, therefore, isAs far as the boundary conditions are concerned, the components of are
4.2. Virtual Variation of the Inertial Work
The virtual variation of the inertial work is where the two dots stand for the second derivative versus time and is the material density. Upon substitution of (13), (28) becomesIn a compact vectorial form, The components of the matrix are where is Kronecker’s delta and
4.3. Governing Equations’ Fundamental Nucleus
The explicit form of the fundamental nucleus of the governing equations is The boundary conditions are For a fixed approximation order, the nucleus has to be expanded versus the indexes and in order to obtain the governing equations and the boundary conditions of the desired model.
5. Closed Form Analytical Solution
The resulting boundary problem is solved via a Naviertype solution. Simply supported beams are, therefore, investigated. The following harmonic displacement field is adopted: where is with representing the halfwave number along the beam axis. is the imaginary unit and the natural frequency. are the maximal amplitudes of the displacement components. The displacement field in (35) satisfies the boundary conditions since Upon substitution of (35) into (33), the fundamental nucleus of the algebraic eigensystem is obtained: The components of the algebraic stiffness matrix are For a fixed approximation order, the eigensystem has to be assembled according to the summation indexes and . Its solution yields as many eigenvalues and eigenvectors (or modes) as the degrees of freedom of the model.
6. Numerical Results and Discussion
Beams made of the aluminium alloy 7075T6 are considered. Mechanical properties are equal to MPa and equal to . The material density is kg/m^{3}. Analyses are carried out considering box and I and Cshaped thinwalled cross sections. The ratio between a representative dimension of the cross section () and the walls’ thickness is . The cross section sides are such that m. A lengthtoside ratio equal to and is considered. Slender and relatively short beams are, therefore, investigated. Classical (bending and torsional) and higher modes up to the axial one are considered. The natural frequency of the axial mode is not included in the investigations since it is accurately obtained by classical EBT and TBT. The natural frequencies are put into the following dimensionless form: As far as validation is concerned, results are compared with threedimensional FEM solutions obtained via the commercial code ANSYS. The quadratic node threedimensional element “Solid186” is used; see ANSYS v10.0 theory manual [31]. The effect of the mesh density is accounted for by considering a fine and a coarse mesh. For the fine mesh, the number of elements along each cross section direction is , whereas elements per side are used in the case of the coarse one. The number of elements along the beam axis is such that the element aspect ratio (defined as the ratio between the axial and each cross section side dimension) is equal to ten. The refined threedimensional solution is addressed as “FEM ” whereas the coarse one is addressed by “FEM .” The modal shapes of the proposed model and the reference threedimensional solutions are compared via visualisation. Mode visualisation was performed within ANSYS postprocessing environment by imposing at each node the displacement components computed by the proposed models through the ANSYS parametric design language command DNSOL; see Madenci and Guven [32]. Since the considered structures are not very complex in terms of geometry, a numerical mode tracking (e.g., by means of the modal assurance criterion, see Allemang [33]) was not necessary. Although the threedimensional FEM solution and the analytical one are different in nature, some considerations about the computational time and effort can be addressed. For the reference FEM simulations, the degrees of freedom (DOFs) over a beam cross section are equal to (fine solution) and (coarse solution) in the case of a box cross section. For a fixed approximation order , the DOFs of the proposed solutions are and they do not depend upon the cross section geometry. In the case of the highest considered expansion order (), they are .
6.1. Box Cross Section
A beam with a box cross section as shown in Figure 2 is considered. Several modes were identified before the axial one. Their number increases as the lengthtoside parameter decreases. Table 2 presents the natural dimensionless frequencies in the case of and equal to one and two. Since the cross section of the beam has two axes of symmetry, two coincident flexural frequencies are present. Their corresponding modes are indicated as Modes I and II. Torsion is addressed as Mode III. Modes I and II for occur before the torsional mode for . For the flexural modes, the hierarchical beam theories match the FEM solution up to four significant digits with an expansion order as low as four. Classical theories yield good results for the flexural modes where a difference is observed only in the fourth significant digit. The natural frequency related to the torsional mode calls for a higher order of expansion. The error for is less than , whereas it is about in the case of a fourthorder model. Consistent with their kinematic hypotheses (the cross section can only translate rigidly on its own plane), classical theories do not yield any torsional mode. The first five dimensionless natural frequencies for the case are presented in Table 3. Modes ordering is referred to as the order of apparition in the model. The same order is observed in the reference FEM threedimensional solution, the modes with being higher than one not considered. Modes I and II are also in this case bending modes. Mode III is a shear deformation on plane as presented in Figure 3. Mode IV is characterised by bending on an opposite direction of two contiguous sides of the cross section; see Figure 4. Mode V corresponds to torsion. Modes III to V call for high expansion in order. Nevertheless, the shear mode differs from the FEM solution by about even for . Rather than further increasing the expansion order, the solution accuracy could be improved by using a local approximation approach in a layerwise sense. This will be matter of future investigations. When the frequencies presented in a table, for a fixed approximation order, are not in an ascending order, the apparition order is not respected. For instance, Mode V is observed before Mode IV for a ninthorder model. Mode swapping is common in loworder solutions; see Giunta et al. [34], and it is due to a higher cross section stiffness. Classical theories are not suitable when because only the two flexural modes are provided. The corresponding values of the frequencies differ from the FEM solution by about four (EBT) and (TBT).
 
a: cross section mesh. b: cross section mesh. c: mode not provided by the theory. 
 
a: cross section mesh. b: cross section mesh. c: mode not provided by the theory. 
(a)
(b)
(a)
(b)
6.2. IShaped Cross Section
Ishaped cross section beams are investigated; see Figure 5. Table 4 presents the first three dimensionless natural frequencies for slender beams. The halfwave number is equal to one and two. The modes are bending on plane , bending on plane , and torsion, respectively. The proposed models all match the FEM solution except for the torsional mode, where the error is about for . The error for the torsion mode considerably increases when the expansion order decreases. It is about , for instance, for a seventhorder model. In the case of , six modes are found before the axial one (considering only equal to one). The first mode is the same as for slender beams. Bending on plane shifts to the third position, the second mode being a torsional deformation. Modes IV to VI present local deformations (“shelllike” behaviour) of the horizontal and vertical cross section elements. These modes are shown in Figures 6–8. The corresponding frequencies are reported in Table 5. The first three frequencies are accurately predicted. The error for the bending mode on plane is about five and in the case of TBT and EBT, respectively. predicts a fourth frequency with a fair accuracy, the error being about . Modes V and VI show the limit of the proposed models and, as a possible solution, a layerwise local approach towards modelling should be used.
