TIMOSHENKO BEAM THEORY EXACT SOLUTION Article Timoshenko Beam Theory: Exact Solution for First-Order Analysis, Second-Order Analysis, and Stability

This paper presents an exact solution to the Timoshenko beam theory (TBT) for first-order analysis, secondorder analysis, and stability. The TBT covers cases associated with small deflections based on shear deformation considerations, whereas the Euler–Bernoulli beam theory (EBBT) neglects shear deformations. Thus, the Euler– Bernoulli beam is a special case of the Timoshenko beam. The momentcurvature relationship is one of the governing equations of the EBBT, and closed-form expressions of efforts and deformations are available in the literature. However, neither an equivalent to the momentcurvature relationship of EBBT nor closed-form expressions of efforts and deformations can be found in the TBT. In this paper, a momentshear forcecurvature relationship, the equivalent in TBT of the momentcurvature relationship of EBBT, was presented. Based on this relationship, first-order and second-order analyses were conducted, and closed-form expressions of efforts and deformations were derived for various load cases. Furthermore, beam stability was analyzed and buckling loads were calculated. Finally, first-order and second-order element stiffness matrices were determined.

beam becomes rigid in shear and the TBT converges towards EBBT. Thus, the Euler-Bernoulli beam is a special case of the Timoshenko beam. First-order analysis of the Timoshenko beam is routinely performed; the principle of virtual work yields accurate results and is easy to apply. Unfortunately, second-order analysis of the Timoshenko beam cannot be modeled with the principle of virtual work. Various studies have focused on the analysis of Timoshenko beams, most of which using numerical methods. Sang-Ho et al. (2019) presented a nonlinear finite element analysis formulation for shear in reinforced concrete Beams; that formulation utilizes an equilibrium-driven shear stress function. Abbas and Mohammad (2013) suggested a two-node finite element for analyzing the stability and free vibration of Timoshenko beam; interpolation functions for displacement field and beam rotation were exactly calculated by employing total beam energy and its stationing to shear strain. Hayrullah and Mustapha (2017) performed a buckling analysis of a nano sized beam by using Timoshenko beam theory and Eringen's nonlocal elasticity theory; the vertical displacement function and the rotation function are chosen as Fourier series. Chen et al. (2018) used the variational iteration method to analyze the flexural vibration of rotating Timoshenko beams; accurate natural frequencies and mode shapes under various rotation speeds and rotary inertia were obtained. Onyia and Rowland-Lato (2018) presented a finite element formulation for the determination of the critical buckling load of unified beam element that is free from shear locking using the energy method; the proposed technique provides a unified approach for the stability analysis of beams with any end conditions. Pavlovic and Pavlović (2018) investigated the dynamic stability problem of a Timoshenko beam supported by a generalized Pasternak-type viscoelastic foundation subjected to compressive axial loading, where rotary inertia is neglected; the direct Liapunov method was used. In stability analysis Timoshenko and Gere (1961) proposed formulas to account for shear stiffness by means of calculation of buckling loads of the associated Euler-Bernoulli beams. Hu et al. (2019) used matrix structural analysis to derive a closed-form solution of the secondorder element stiffness matrix; the buckling loads of single-span beams were also determined.
The momentcurvature relationship, one of the governing equations of the EBBT, has no equivalent in the TBT literature. Furthermore, closed-form expressions of efforts and deformations are not common. In this paper a momentshear forcecurvature relationship (MSFCR), the equivalent in TBT of the momentcurvature relationship of EBBT, was presented. The relationship between the curvature, the bending moment, the bending stiffness, the shear force, and the shear stiffness was then described. Based on MSFCR, closed-form expressions of efforts and Substituting Equation (3) into Equation (6) yields the following: Equation (7) yields the following momentshear forcecurvature relationship (MSFCR) combining bending and shear: For a uniform beam along segments, substituting Equation (1) into Equation (8) yields, For a tapered beam, substituting Equation (1) into Equation (8) yields, In the case of non-uniform heating, Equation (9) becomes: where  T is the coefficient of thermal expansion, T = T bb -T tb is the difference between the temperature changes at the beam's bottom fibers (T bb ) and top fibers (T tb ), and d is the height of the beam.
Combining Equations (1) and (5) yields the rotation angle as follows: Equations (5) to (12) apply as well in first-order analysis and in second-order analysis.  The bending shear factor is defined as follows:

Summary of Timoshenko and Euler-Bernoulli beam equations
(13)

Governing equations
The application of Equation (2) yields the following formulation of the bending moment: Substituting Equation (2) into the MSFCR (Equation (9)) for a uniform beam and integrating twice yields: The shear forces and rotation angles are determined using Equations (1) and (14), and Equations (12), (14), and (15), respectively. (15a) The integration constants C i (i = 1, 2, 3, and 4) are determined using the boundary conditions and continuity equations and combining the deflections, the rotation angles, the bending moments, and the shear forces.
The element stiffness matrix in local coordinates of the Timoshenko beam is denoted by K Tbl . The relationship between the aforementioned vectors is as follows: Applying Equation (15c) with the distributed load q = 0, yields the following: Considering the sign conventions adopted for bending moments and shear forces in general (see Figure 1) and for bending moments and shear forces in the element stiffness matrix (see Figure 2), we can set following static compatibility boundary conditions using Equation (19): Considering the sign conventions adopted for the displacements and rotations in general (see Figure 1) and for displacements and rotations in the member stiffness matrix (see Figure 2), we can set following geometric compatibility boundary conditions using Equations (20)    The vectors of Equations (16) and (17) become The element stiffness matrix becomes: where (34)

Beam resting on an elastic Winkler foundation
For a beam resting on a Winkler foundation with stiffness K w , Equation (2) of static equilibrium becomes, Differentiating Equation (35) twice with respect to x, combined with the MSFCR (Equation (9)), yields: The solution of Equation (36) yields the formulation of M(x) with four integration constants (A 1 , B 1 , C 1 , D 1 ). The shape of M(x) depends on the parameter k w -4/ 2 as follows: for k w -4/ 2 > 0, k w -4/ 2 = 0, and k w -4/ 2 < 0, respectively . Combining M(x) with Equation (35) yields the following equations: The expression of M(x) being known, applying Equation (1) yields the shear force, and the combination of M(x) with Equations (12) and (38) yields the rotations of the cross section.

Second-order analysis of uniform Timoshenko beams
The analysis conducted in this section holds for compressive forces smaller (absolute values) than the buckling load.

Governing equations
A uniform beam is analyzed, an elastic Winkler foundation not being considered. The axial force N (positive in tension) is assumed to be constant in segments of the beam. The equation of static equilibrium is as follows: Combining Equation (39) with the MSFCR (Equation (9)) yields the following: The solution of Equation (40) yields M(x), which contains the integration constants C 1 and C 2 . The combination of M(x) with Equation (39) yields following equations: The transverse force T(x) is determined as follows: Substituting Equations (1) and (41) into Equation (43) yields the transverse force T(x) as follows:  (41)   Combining Equations (39) and (49) yields the following: The combination of Equation (12) for the rotation of the cross section with Equations (49) and (51) yields the following, with the parameter 2 defined as shown: By applying Equations (1), (43), (49), and (51), the transverse force T(x) yields: Considering the static compatibility boundary conditions (Equations (21) to (34)), whereby the shear forces are replaced by the transverse forces, the geometric compatibility boundary conditions, and Equations (49) to (55), the following equations are obtained: The combination of Equations (18), (45), (46), and (56) yields the element stiffness matrix as follows: (57) 1 1 1 1 1 1 1 The development of Equation (57) yields the following formulation: Accordingly, the element stiffness matrix of the beam having a hinge is as follows: (59) Case 2: Tensile force or compressive force with k  -1 The solution of Equation (48) is as follows, with the parameter 3 defined as shown: Combining Equations (39) and (60) and integrating with respect to x yields the following: Combination of Equation (12) for the rotation of the cross section, with Equations (60) and (62) yields the following: The parameter  4 (Equation (65)) has a positive value in the case of tension and a negative value in the case of compression with k  -1. Applying Equations (1), (43), (60), and (62), the transverse force T(x) yields the following: K 11 = -k 2 sin 1 / 1   11  12  11  12  3  2  3  2   22  12  24  2   11  12  3  2 22 Considering the static compatibility boundary conditions (Equations (21) to (24)), whereby the shear forces are replaced by the transverse forces, and the geometric compatibility boundary conditions, and Equations (60) to (66), the element stiffness matrix can be expressed as follows: The development of Equation (67) yields the following formulation: Accordingly, the element stiffness matrix of the beam having a hinge is as follows:

Analysis of some load cases
Uniformly distributed load p.

Beam resting on an elastic Winkler foundation with an axial load
In the equation of static equilibrium, the axial force N (the value of which is positive in tension) is assumed to be constant in segments of the beam. The stiffness of the Winkler foundation is denoted by K w .
The combination of Equation (72) with the MSFCR (Equation (9)) yields the following: Differentiating twice both sides of Equation (73) with respect to x and combining the result with the MSFCR (Equation (9)) yields the following: The solution of Equation (74)

yields the formulation of M(x) with four integration constants. From M(x), combined with
Equations (1) and (73), the shear force V(x) and the deflection w(x) can be deduced. The application of Equations (12) and (43) yields the transverse forces T(x) and the rotations of the cross section (x).

Stability of the Timoshenko beam
The formulations of bending moment M(x), deflection w(x), rotation of the cross section (x), and transverse force T(x) are used to determine the buckling load.
For k  -1, Equations (49) to (55) are considered to satisfy the boundary conditions and continuity conditions.
For k  -1, Equations (60) to (66) are considered to satisfy the boundary conditions and continuity conditions.
The resulting eigenvalue problem is solved to determine the buckling loads.

TIMOSHENKO BEAM THEORY EXACT SOLUTION
The buckling load N cr is defined as follows: (75) Values of the buckling factor  are listed in Table 4. Closed-form expressions of the matrices expressing the boundary conditions are presented in Appendix D. To determine the buckling loads the determinants of those matrices are set to zero. Closed-form expressions of the buckling factors for a pinned-pinned beam, a fixed-free beam, and a fixed-fixed beam are also presented in Appendix D. We recall that the exact values of the buckling factors  for  = 0.0 (corresponding to the Euler-Bernoulli beam) for the support conditions SS-SS, F-SS, F-FR and F-F are 1.00, 0.700, 2.00, and 0.500, respectively.

Beams subjected to uniformly distributed load and an axial force
Let us calculate the responses of a uniform pinnedpinned beam subjected to a uniformly distributed load and an axial force, as shown in Figure 7.
The results of this study are in good agreement with those of Hu et al. (2019).
In fact, both formulas are identical since following equivalences exist between the parameters considered by Hu et al.
The equations for the determination of the buckling loads (see Appendix D) are also identical to those of Hu et al.

Conclusion
The momentshear forcecurvature relationship presented in this study enabled the derivation of closed-form solutions for first-order analysis, second-order analysis, and stability of Timoshenko beams. The results showed that the calculations conducted as described in this paper are exact. Closed-form expressions of efforts and deformations for various load cases were presented, as well as closed-form expressions of second-order element stiffness matrices (the axial force being tensile or compressive) in local coordinates. The determination of element stiffness matrices (ESM) enables the analysis of systems with the direct stiffness method. We showed that ESM can also be determined by the presence of hinges.

Calculation of bending moments and shear forces:
Influence of tensile force: With increasing tensile force, bending moments decrease (in absolute values), and with increasing bending shear factor, bending moments decrease (in absolute values).
Influence of compressive force: With increasing compressive force, bending moments increase (in absolute values), and with increasing bending shear factor, bending moments also increase (in absolute values).
Stability of the beam: With increasing bending shear factor, the buckling load decreases.
The following aspects not addressed in this study could be examined in future research:  Analysis of linear structures, such as frames, through the transformation of element stiffness matrices from local coordinates in global coordinates.
 Second-order analysis of frames free to sidesway with consideration of P- effect.
 Use of the direct stiffness method, since element stiffness matrices are presented.
 Analysis of positions of discontinuity (interior supports, springs, hinges, abrupt change of section), since closedform expressions of bending moments, shear or transverse forces, rotation angles, and deflection are known.


Beams resting on Pasternak foundations, the Pasternak soil parameter can be considered as a tensile force.

Supplementary Materials:
The following files were uploaded during submission:  "Deflection calculation of a pinnedpinned beam using the principle of virtual work,"  "Analysis of a fixedpinned beam under concentrated load,"  "Buckling analysis of a pinnedpinned beam,"  "Buckling analysis of a fixedpinned beam,"  "Buckling analysis of a fixedfree beam,"  "Buckling analysis of a fixedfixed beam,"  "Bending analysis of a pinnedpinned beam subjected to a uniformly distributed load and an axial force,"  "Bending analysis of a fixedpinned beam subjected to a concentrated load and an axial force," and  "Second-order element stiffness matrix."

Conflicts of Interest:
The author declares no conflict of interest.

Appendix C Fixedpinned Timoshenko beam subjected to a concentrated load
Equations (19), (20) For calculation with the principle of virtual work the fixed-end moment and the moment under the load can be expressed as follows: (C3) For other support conditions (pinnedpinned, fixedfree, and fixedfixed) the boundary conditions (Equations (C2) a , (C2) b , (C2) g , and (C2) h ) are modified accordingly, and the bending moment and deflection curves are deduced using Equation (15c) with q = 0.

Appendix D Buckling loads of single-span beams with various support conditions
Compressive force with k  -1: The eigenvalue problem led to following equations: sinh 3 = 0, sinh 3 - 4 cosh 3 = 0, cosh 3 = 0, and 2 -2cosh 3 +  4 sinh 3 = 0 for a pinned-pinned beam, a fixedpinned beam, a fixed-free beam, and a fixed-fixed beam, respectively. Recalling that  4 is negative, those equations have no nontrivial solutions; consequently, no buckling occurred due to the action of a compressive force with k  -1.

Compressive force with k > -1
Pinned   with the matrix and the equations as follows: (D4)

Appendix E Timoshenko beams subjected to various load cases and an axial force
The compressive forces are assumed smaller (absolute values) than the buckling load.
Case 1: Compressive force with k  -1. The bending moment and the deflection curve are denoted by M c1 (x) and w c1 (x), respectively. Let us recall the parameters  1 and  2 defined in Equations (50) and (54), respectively.

(E1)
Case 2: Tensile force or compressive force with k  -1. The bending moment and the deflection curve are denoted by M c2 (x) and w c2 (x), respectively. Let us recall the parameters  3 and  4 defined in Equations (61) and (65), respectively. (E2) The parameter  4 (Equation (65)) has a positive value in the case of tension and a negative value in the case of compression with k  -1.