SciELO - Scientific Electronic Library Online

vol.11 issue12Study of reflection and transmission of plane waves at thermoelastic-diffusive solid/liquid interfaceParametric study on average stress-average strain curve of composite stiffened plates using progressive failure method author indexsubject indexarticles search
Home Pagealphabetic serial listing  

Services on Demand




Related links


Latin American Journal of Solids and Structures

On-line version ISSN 1679-7825

Lat. Am. j. solids struct. vol.11 no.12 Rio de Janeiro  2014 

General exact harmonic analysis of in-plane timoshenko beam structures



C. A. N. Dias

Group of Solid Mechanics and Structural Impact Department of Mechatronics and Mechanical System Engineering University of São Paulo - São Paulo, Brazil E-mail:




The exact solution for the problem of damped, steady state response, of in-plane Timoshenko frames subjected to harmonically time varying external forces is here described. The solution is obtained by using the classical dynamic stiffness matrix (DSM), which is non-linear and transcendental in respect to the excitation frequency, and by performing the harmonic analysis using the Laplace transform. As an original contribution, the partial differential coupled governing equations, combining displacements and forces, are directly subjected to Laplace transforms, leading to the member DSM and to the equivalent load vector formulations. Additionally, the members may have rigid bodies attached at any of their ends where, optionally, internal forces can be released. The member matrices are then used to establish the global matrices that represent the dynamic equilibrium of the overall framed structure, preserving close similarity to the finite element method. Several application examples prove the certainty of the proposed method by comparing the model results with the ones available in the literature or with finite element analyses.

Keywords: exact harmonic analysis; Laplace transform; Timoshenko beam; dynamic stiffness matrix; rigid offsets; end release.




j imaginary number

t time variable

χ flexible span axial coordinate

ω excitation circular frequency

y flexible span transverse


Ε elastic modulus

ρ mass density

Α cross section area

Ι bending inertia

m mass per unit length

G shear modulus

υ Poisson's ratio

Αs shear area

Ρ static axial load

r rotatory inertia per unit


cΙ internal damping

αΑ slope of the distributed axial load

cΕ external damping

αL slope of distributed

Transverse load

bΑ distributed axial load χ = 0 at χ = 0

LΤ total member length

L member flexible span length

ΡΑΙ, ΡΑJ distributed axial load at node I and at node J, respectively

ΡLI, ΡLI dist ributed lateral load at node I and at node J, respectively

(αI, bI), (αJ, bJ)dimensions of the rigid offset attached at node I and at node J, respectively

LI,LJ length of the rigid offset attached at node I and at node J, respectively

MI, MJ mass of the rigid offset attached at node I and at node J, respectively

MRI, MRJ rotatory inertia of the rigid offset attached at node I and at node J,




υ (x,t) axial displacement

VB (x,t) bending deflection

V (x,t) total deflection

VS (x,t) shear deflection

θ (x,t) bending


normal stress

normal force

bending stress

shear force

shear stress



amplitude of any of the above functions

Laplace transform of



end displacement and end force vectors of a member flexible span, respectively

end displacement and end force vectors at member nodes, respectively

equivalent load vectors at member flexible span and at member nodes, respectively

member dynamic stiffness, rotation, and connection matrices, respectively

global vectors of nodal displacements and of nodal forces, respectively

global vector of the equivalent loads and global dynamic stiffness matrix, respectively

global matrices of the nodal concentrated masses and springs, respectively



Many modern structures are formed by beam elements. These skeleton like structures are subjected to static and dynamic loads. Their beam members can be of various sizes, including beams with small length to beam height ratio. For the analysis of these structures, it is important to use a more refined beam theory, where the assumption of the cross section to remain plane is not enforced. Besides, harmonic loads can be of high frequency, when then it is important to keep in the beam model the cooperation of the rotatory inertia to the overall structure response. This motivates the use of the Timoshenko beam theory to obtain the damped steady state response for general plane frames subjected to harmonically time varying external forces. The study reported here concerns with an exact harmonic analysis using the classical Dynamic Stiffness Matrix, DSM. The problem at hand is non-linear and transcendental with respect to the excitation frequency [Howson and Williams (1973)]. The approach used to solve it is by the use of the Laplace transform.

Focusing on the calculation of natural frequencies, Howson and Williams (2003) present a formulation based on the classical DSM obtained by a set of decoupled fourth order partial differential equations (PDE) for the unknowns deflection and rotation of the beam cross section. Dias and Alves (2009) also derive the DSM via the same procedure but reaches an improved formulation, argued to be more suitable for the eigenproblem solution. It has been noticed [Schanz and Antes (2002)] that the dynamic analyses of beams can be performed by decoupling deflection and rotation. In the study described here, the original coupled PDEs, combining deflection, rotation, bending moment and shear force, are directly employed in order to reach the DSM formulation as well as the formulation for the equivalent load vector that arises due to distributed loads on the beams and frames.

We remark that, to be considered exact, a solution must adopt no assumptions other than those of the beam theory itself. Hence, if the mode superposition method is employed, the result is not exact once it is affected by series truncation error. For this reason, many good publications dedicated to dynamic frame response calculation do not conform to the requirements of an exact solution and therefore are not considered here. Concerning the dynamic response of frames, Table 1 compares different solution methods found in the literature [Abu-Hilal (2003), Foda and Albassam (2006), Gürgöze and Erol (2001), Lin (2008), Loudini et. al. (2006), Majukt (2009), Mehri et. al. (2009), Saeid (2011), Schanz and Antes (2002), Tang (2008)]. In the present context, it is imperative that the PDE representing the governing equations of a given in-plane structure subjected to harmonic forces has to be solved exactly. This can be achieved by using Laplace transform [Abu-Hilal (2003), Foda and Albassam (2006), Loudini et. al. (2006), Saeid (2011)] and/or Green functions [Abu-Hilal (2003), Tang (2008), Davar and Rahmani (2009)]. Gürgöze and Erol (2001) use a distinct method based on a receptance matrix but cannot be considered exact since series solution is employed.

As seen in Table 1, only the solutions given in (Shanz and Antes, 2002; Foda and Albassam, 2006; Lin, 2008; Majkut, 2009; Saeid, 2001) could be generalized in order to solve models of arbitrary numbers of beams and boundary conditions. If further features are to be considered, e.g.

concentrated masses and springs, these solutions would be limited to the ones developed in (Foda and Albassam, 2006; Majkut, 2009). Even so, these references are dedicated to solve single and/or in-line beams. When considering the capability to solve framed structures with the features of concentrated masses and springs and rigid bodies attached to the members ends, only Seeid (2001) can be highlighted. These features are not taken into account by Antes et al. (2004), which also deals with harmonic loads applied to Timoshenko frames. Mei (2008, 2012) present an interesting model that considers in-plane vibration of some restricted forms of frames. The analyses presented in these references are developed, as here, along the Timoshenko bending theory. In contrast, the present paper is more general inasmuch as it handles any shape of portal planar frame subjected to harmonic loads and it solves the equations of motion via the Laplace transform.

Except for the case of transient analysis, which is beyond the scope of this study, the solution of the present investigation is unique in the sense that attends to all of the requisites listed in

Table 1. Indeed, to best of the author knowledge, no other publication in the context of harmonic analysis pays attention to effects like rigid offset and end release. Few publications (Abu-Hilal, 2003; Tang, 2008; Saeid, 2001) have considered distributed loads and even fewer have incorporated damping (Loudini, 2006; Abu-Hilal, 2003) in their analysis.

The framed structure under consideration is composed by in-plane members, which are connected to each other through previously defined nodes. Each member has a flexible span, which obeys the Timoshenko beam theory, and may also have rigid bodies attached to its ends,

where forces can be released. Fig. 1 depicts the internal forces and displacements for the flexible span, while Fig. 2 shows the end forces, end displacements, distributed loads, and attached rigid bodies for a typical member.



Since the rigid bodies may or may not have mass, herein they are alternatively named rigid offsets. Taking into account external and internal damping effects, Section 2 is dedicated to presenting the basic equilibrium equations of a flexible span, in both axial and transverse directions. For the transverse equilibrium, the Timoshenko beam theory (Howson and Williams, 1973; Dias and Alves, 2009; Schanz and Antes, 2002; Loudini et al., 2006) is applied in order to consider the effects of shear deflection and distributed rotatory inertia. In Section 2.2, under the assumption in Section 2.1 that all the excitations vary harmonically with the same frequency, the steady state

response is formulated in terms of displacement and force amplitudes, which are then subjected to the Laplace transform. Although this assumption reduces the application scope, it characterizes the well-known harmonic analysis, useful to identify harmful vibrations due to resonance occurrences. On the other hand, considering that a given beam like structure is linear, the proposed method can be extended to a more general excitation type, with more than one distinct excitation frequency.

The first order governing equations of Section 2 are a coupled PDE system that might be manipulated in order to obtain uncoupled PDEs for the displacements (Howson and Williams, 1973; Foda and Albassam, 2006). However, this operation would produce fourth order PDEs that are less suitable for Laplace transform (Schanz and Antes, 2002). Hence, in Section 3, the equations of Section 2 are directly subjected to Laplace transforms such that, after some easier but laborious algebraic operations, decoupled results for displacements and forces could be obtained in the Laplace domain.

After applying inverse Laplace transforms to the results of Section 3, Section 4 gives the desirable formulations for the internal displacements and forces of the member flexible span, which are expanded to the model global matrices in Section 5. As it is shown there, from these formulations it is possible to obtain the relation between end forces and end displacements such that the member dynamic stiffness matrix and the member equivalent load vector can be defined. Then, using a similar technique to the finite element method, it is shown how these member matrices are employed to establish the global matrices that represent the dynamic equilibrium of the framed structure.

The solution of this equilibrium system directly gives the nodal displacements, from which easily follows member stresses. Some careful chosen examples in Section 6 illustrate the major features of the present method.



According to the Timoshenko beam theory (Howson and Williams, 1973; Dias and Alves, 2009;

Schanz and Antes, 2002; Loudini et al., 2006) and considering the internal forces and displacements in Fig. 1, with the distributed linear varying external forces in Fig. 2.b, the following equations are obtained, with the prime denoting the derivative wrt position, χ , overdo denoting derivative wrt time, t , and all the loads varying in time with the same circular frequency ω :

Normal stress, σN



Normal force, FN



Axial equilibrium



Mean shear stress, τ



Shear force, FS



Transverse force equilibrium



Bending stress, σB



Bending moment, MB



Moment equilibrium



Here, E and G are the elastic and shear moduli, u and v are the axial and total transverse beam displacement, vs is the transverse displacement due to shear, A and As are the total cross-section and shear areas, αAand αL are the slope of the distributed axial and transverse load, bA and bL are the distributed axial and transverse load, m is the mass per unit length, ρ and P are distributed axial loads.

In these equations, the internal and external damping coefficients, cI and cE, respectively, are defined in order to reproduce the Rayleigh formulation, where cErepresents the mass proportional coefficient and cI the stiffness proportional coefficient. From this, by using the orthogonal condition of the natural modes, it is possible to establish the following expression for the modal damping ratio ζI :



Meaningful values for cI and cE can be obtained whenever a pair of modal damping ratio values is known. Therefore, after applying Eq. (10), these damping coefficients can be computed by:



where the natural frequencies ( ω1, ω2 ) and the damping ratios ( ζ;1, ζ;2 ), for the lowest first and second modes, are supposedly known. On the other hand, by imposing that both Rayleigh damping coefficients must be non-negative and that the peak frequencies ωpi given by:



must be real, the following restrictions must be obeyed:



The external damping coefficient cE accounts for the environment viscous actions expressed by the forces (x,t) and cv(x,t) , so that:



The internal damping coefficient cIaccounts for the energy dissipation due to the structure deformation rate and can be related with the Kelvin-Voight damping coefficient KV (Shanz and Antes, 2002):



A more complete characterization of damping in dynamic structural systems can be found in Oliveto, and Greco (2002), where it is shown that the Rayleigh coefficients are independent of the boundary conditions, no matter whether the system is continuous or discrete.

2.1 Steady state response

For the steady state response, all the previous defined internal forces and displacements can be written according to the general equation



with time and space being decoupled.

Performing time derivatives in equations (1-9), the amplitudes are related by:



where, for instance, , is the amplitude of the beam total deflection.

2.2 Laplace Transform

Applying the Laplace transform over the previous defined amplitudes:



the differential equations (16) are substituted by the following simple algebraic expressions:



It is adopted in the sequence, according to Fig. 2.a, the convention for the end forces and end displacements as



Now, solving Eqs. (18.a) and (18.b) for the axial displacement results in:









Total deflection can be obtained after considerable algebraic manipulation of Eqs. (18.c-f), yielding:






when adopting






Analogously, solving Eqs. 16(c-f) for the shear force and bending rotation results in:






The Laplace transform of the bending moment can be obtained by inserting Eq. (28) into Eq. (18.e). The bending moment amplitude,, comes from Eq. (16.e) once the rotation amplitude, , has been obtained.

A common practice in the specialized literature is to decouple Eqs. (16) in order to obtain a set of differential equations involving only one unknown function. For the Timoshenko theory, complicated terms involving fourth order derivatives will appear, so the use of the Laplace transform is troublesome. To avoid this, here the Laplace transform was directly applied to the coupled first order Eqs. (16) by applying the transform to second order coupled equation

involving only deflection and rotation. As a remark, Schanz and Antes (2002) also avoids to work with fourth order derivatives.



By applying the inverse transform to Eq. (20) for the axial displacement and then using the result in Eq. (16.a) for the normal force, these corresponding amplitudes are then given, respectively, by:



where, knowing that






it holds that

1. For n = 0,2,4,6.....



2. For n = 1,3,4,5.....



3. For n = -2,-4,-6 .....



4. For n = -1, -3,-4, ....



As for the total deflection, it can be obtained by applying the inverse transform to Eq. (23), yielding












The rotation due to bending is obtained by substituting Eq. (23) into Eq. (28) and applying the inverse transform so that the amplitude is:









Finally, the shear force amplitudes come from the inverse transform of Eq. (27), resulting in



while the bending moment amplitude comes from the use of Eq. (37) in expression (16.e),




4.1 Flexible span equilibrium

Applying the following boundary conditions for the end displacement (see Fig. 2.a):



after laborious algebraic work, the equations of the previous section can be rearranged in the following expressions for the internal displacements:



with the functions on the right hand sides (from being defined in Appendix A. Analogously, applying the following boundary conditions for the end forces:



the internal forces are given by:



with the functions on the right hand sides being defined in Appendix A. Now, defining the vectors:






and using the previous equations (43.a) for the axial displacement and (45.a) for the normal force, the axial equilibrium is expressed by:






Equation (48), after multiplication by the inverse of φPA , gives the following classical matrix equilibrium equation:



with the equivalent lateral load vector and the dynamic stiffness matrix given by, respectively,






Now, combining both axial and lateral matrix equilibrium equations, it follows that:



with the end displacements and end forces vectors being given, respectively, by



The equivalent load vector is



and the dynamic stiffness matrix is



4.2 Rigid offset

Observing Fig. 2, the conversions from the flexible span ends to the member ends are given by:






Now, collecting from expressions (64, 65) and substituting into Eq. (59), member equilibrium requires that:



where, for the member equivalent load vector and the member dynamic stiffness matrix, respectively, we have:



4.3 End force release

Let IP denote a permutation matrix that moves, to the last line, the degree of freedom to be released. If, for any end displacement vector Q, it is imposed that the corresponding end force must be null, it can be written that:



By doing so, the released versions for the member equivalent load vector and for the member dynamic stiffness matrix are, respectively:



with the partitions coming from:



Finally, the displacement of the released degree of freedom can be obtained from:




The member equilibrium equations seen in the previous sections can be also expressed as



with the subscript i added to identify the i-th member of the model and the dependency on the excitation frequency ω is explicitly indicated. Then, the usual FEM model assembly leads to the global dynamic equilibrium written as



where is the vector of the applied nodal loads, is the vector of the nodal displacements, both in global coordinates, and



is the global equivalent load vector, withdenoting the member rotation matrix andenoting the connection matrix that relates the member degrees of freedom with the model degrees of freedom. The global dynamic stiffness matrix is expressed by:



with denoting the matrix of the nodal concentrated springs and denoting the matrix of the nodal concentrated masses. Once the global equilibrium equation (80) has been assembled, the solution is obtained by:



with the backslash operator denoting an appropriate solution method, based, for instance, on the Gauss elimination. Once the vector has been obtained, the member end displacements are evaluated by:



and the member end forces by:




The exact harmonic analysis proposed herein was implemented in the standalone application VIHAND, developed in the MATLAB platform. This code is part of the VIANDI computational system for the exact static and dynamic analysis of skeleton-like structures (Dias, 2011) and freely available at In order to explore the above exact solution, some examples were judiciously chosen and compared with other analytical or numerical approaches from the literature, as now described.

6.1 Euler-Bernoulli beam

When shear deflection and rotatory inertia are disregarded, the Timoshenko theory must reproduce the Euler-Bernoulli theory. Therefore, the aim of this section is to verify the present method for the particular case when by comparing it with results from other authors. The examples given in Tables 2 and 3, where no damping is presented, show a very good agreement with the exact solution given in (Abu-Hilal, 2003), which uses Green functions obtained by the Laplace transform.

In both tables, additional results obtained by the approximate mode superposition method for different number of modes are also listed. Clearly, it can be seen that the mode superposition series solution improves with increasing number of modes and that the convergence is quite good for the deflection but not so good for the shear force. This confirms the well-known fact that the mode superposition method lacks accuracy for predicting variables involving derivatives wrt position.

For the damped cases in Table 4, Abu-Hilal (2003) does not supply written tabular results. For this reason, the results were obtained by using the formulae presented there. Additionally, a FEM model was employed to solve the problem. 1000 beam elements have to be used to model the beam so to reach the same accuracy as that given by the present method. As can be observed in the table, the agreement of the three methods is quite good.

Based on a receptance matrix, Gürgöze and Erol (2001) offer a series solution for the example in Table 5 as well as an exact solution based on boundary value formulation. Again, the agreement with the present method is very good.

6.2 Timoshenko beam

For the examples with no damping listed in Table 6, Foda and Albassam (2006) use Green functions to provide exact solutions for Timoshenko beams with attached masses and springs. Unfortunately, the authors do not supply written tabular results, so only approximations could be collected from the graphics presented there. As seen in the table, the reliability of the present

method could be confirmed by comparison with a FEM model with 1000 beam elements. By adding damping effects, the same problem has the frequency response depicted in Fig. 3, where a very good agreement with the FEM solution can be observed.

6.3 Elastic robot manipulator

The authors did not find in the open literature examples of the exact behavior of more complex skeleton structures. They would allow exploring further the resources in VIHAND and the theory here presented. However, a more refined model representing a robot manipulator is suggested in Fig. 4, whose free vibration response is studied in (Dias, 2010). Here, the effects of rigid offset, bending release, skewed edge, static and dynamic loads, concentrated mass, long and short beams and an external spring, besides damping, are all present. The results obtained by the theory in this article are then compared with the FEM, which is the only alternative that permits to solve such a problem.

It is interesting to notice that VIHAND is capable of dealing with the interaction between static and dynamic loads. In the case of this robot structure, when the horizontal arm of the manipulator receives a static compression load, the natural frequencies of the structure decreases. This leads to a harmful vibration due to the resonance between the actuator excitation frequency and the first natural mode. As can be seen in Fig. 5, between 12 and 13 Hz the dynamic stress increases more than twofold, a fact that can severely compromise the structure fatigue life.

The results of the VIHAND and ANSYS programs agree very well. However, it must be pointed out that the ANSYS model needed more than 240 beam elements to achieve the same accuracy as that obtained by the exact solution of the program VIHAND, which employs only four elements. Two aspects of the present model are worth considering. First, the method is not suitable for spatial frames analyses, at least in an exact way. The reason is that, for these structures, it is likely that torsion-flexure coupling will be important and torsion is not here taken into account.

The consideration of this coupling is far from simple and no attempt was made in the present study to handle it. The second point to highlight is that this model does not intend to be an alternative to the finite element method. But, being possible to obtain an exact solution, it is certain that finite element procedures can benefit from this study inasmuch as the solution procedure here outlined can serve as a reference for this and other numerical procedures.



By applying a Laplace transform to the coupled first order PDEs, which combine displacement and force amplitudes, this work consolidates an exact solution method for the frequency response of in-plane framed structures. The method, based on the Timoshenko beam theory, is quite general since it allows the use of concentrated masses and springs, rigid bodies, end force release, skewed edges, internal damping and external damping as well as concentrated and distributed harmonic loads.

The reliability of the proposed method could be tested by comparing it against simple beam examples found in the literature. The method was also fully verified by means of a more complex model involving the majority of its capabilities and comparing the results with the FEM. Although a priori developed for a single excitation frequency, the method can be easily extended for more complex and arbitrary time varying excitations by using Fourier series and/or Fourier transform. Hence, the method is useful for spectral analysis implementation in the context of random vibrations.



Abu-Hilal M (2003) Forced vibration of Euler-Bernoulli beams by means of dynamic Green functions, Journal of Sound and Vibration, 267: 191-207.         [ Links ]

Antes H, Schanz M and Alvermann S (2004) Dynamic analyses of plane frames by integral equations for bars and Timoshenko beams, Journal of Sound and Vibration, 276: 807-836.         [ Links ]

Dias CAN (2010) Power Secant Method applied to natural frequency extraction of Timoshenko beam structures, Latin American Journal of Solids and Structures, 7: 307-333.         [ Links ]

Dias CAN (2011) Computational system for dynamic analysis of in-plane beam structures, ISBN 978-85-900395-0-7, published online at         [ Links ]

Foda MA and Albassam BA (2006) Vibration confinement in a general beam structure during harmonic excitations, Journal of Sound and Vibration, 295: 491-517.         [ Links ]

Gürgöze M and Erol H (2001) Determination of frequency response function of a cantilevered beam simply supported in-span, Journal of Sound and Vibration, 247: 372-378.         [ Links ]

Howson WP and Williams FW (1973) Natural frequencies of frames with axially loaded Timoshenko members, Journal of Sound and Vibration, 26: 503-515.         [ Links ]

Lin HY (2008) Dynamic analysis of multi-span uniform beam carrying a number of various concentrated elements, Journal of Sound and Vibration, 309: 262-275.         [ Links ]

Loudini M, Boukhetala D, Tadjine M and Boumehdi MA (2006) Application of Timoshenko beam theory for deriving motion equations of lightweight elastic link robot manipulator, ICGST-ARAS Journal, 5 II: 11-18.         [ Links ]

Majkut L (2009) Free and forced vibrations of Timoshenko beams described by single difference equation, Journal of Theoretical and Applied Mechanics, 47: 193-210.         [ Links ]

Mehri B, Davar A and Rahmani O (2009) Dynamic Green function solution of beams under a moving load with different boundary conditions, Mechanical Engineering Transactions, 16: 273-279.         [ Links ]

Mei C (2008) Wave analysis of in-plane vibrations of H- and T-shaped planar frame structure, ASME Journal of Vibration and Acoustics, 130: 1-10.         [ Links ]

Mei C (2012) Wave Analysis of In-Plane Vibrations of L-Shaped and Portal Planar Frame Structures, ASME Journal of Vibration and Acoustics,134: 1-12.         [ Links ]

Oliveto G and Greco A (2002) Some observations on the characterization of structural damping, Journal of Sound and Vibration, 256: 391-415.         [ Links ]

Saeid AA (2011) Dynamic stiffness matrix and load functions of Timoshenko beam using the transport matrix, Computers & Structures, 79: 1175-1185.         [ Links ]

Schanz M and Antes H (2002) A boundary integral formulation for the dynamic behavior of a Timoshenko beam, Electronic Journal of Boundary Elements, BETEQ 2001, 3: 348-359.         [ Links ]

Tang B (2008) Combined dynamic stiffness matrix and precise time integration method for transient forced vibration response analysis of beams, Journal of Sound and Vibration, 309: 868-876.         [ Links ]



Received 15.05.2014
In revised form 29.06.2014
Accepted 18.08.2014
Available online 26.09.2014



Appendix A

For the beam flexible span in Fig. 2, the following defines the functions needed to calculate the internal displacement and force amplitudes by Eqs. (43-45). Function is defined by Eqs. (33) and by Eqs. (36).

A.1 Axial displacement



A.2 Normal force



A.3 Total deflection



A.4 Bending rotation



A.5 Shear force



A.6 Bending moment


Creative Commons License All the contents of this journal, except where otherwise noted, is licensed under a Creative Commons Attribution License