Levinson Beam Theory: An Analytical Solution to the First-Order Analysis of Uniform Beams with Simply Symmetrical Cross-sections

1. Introduction

The Euler−Bernoulli Beam Theory (EBT) also called “Classical Beam Theory” is the first beam theory formulated for the analysis of beams. It was first enunciated circa 1750 and did not consider shear deformations. Therefore the deflections predicted were only accurate for slender beams: in case of short or thick beams the deflections were underestimated. To overcome the drawbacks posed by EBT a first-order shear deformation theory was proposed 1921 by Timoshenko [1]. Shear deformations were considered but shear coefficients had to be defined. Moreover, the Timoshenko beam theory predicted constant shear stresses in the cross-section, this violating the shear stress free surface conditions. To deal with the problems of the first-order shear deformation theory, higher-order shear deformation theories were developed. In this respect Levinson [2] proposed 1981 a third-order shear deformation theory for beams with rectangular cross-sections whereby the equations were derived using a vectorial formulation. The Levinson beam theory (LBT) did not use shear coefficients and the shear stress free surface condition was satisfied. The Reddy beam theory [3] proposed in 1984 used the same displacement field as LBT but the governing differential equations were based on variationnaly derived equations of equilibrium. Other research papers on shear deformation theories can be mentioned. Golbakhshi et al. [4] presented a revised couple stress model similar to the classical Timoshenko beam theory according to the third-order Reddy–Levinson beam theory. Faghidian et al. [5] conceived and applied a stationary variational framework of the Timoshenko–Ehrenfest beam, founded on the elasticity theory, to set forth the variationally consistent shear coefficient of a prismatic beam of a solid rectangular cross-section. Chen et al. [6] presented accurate values for the dynamic stiffness matrix coefficients of Levinson beams under axial loading embedded in a Winkler–Pasternak elastic foundation, while Ike [7] presented the determination of natural transverse vibration frequencies of Euler-Bernoulli beams on Winkler foundations using the Generalized Integral Transform Method. Schoeftner [8] proposed an analytical refined beam model taking into account the axial, the shear and the transverse normal stress distributions whereas Karttunen et al. [9] provided a consistent variational formulation for the static Levinson beam theory. The variational formulation was carried out by employing the principle of virtual displacements whereby the external virtual work done by the stresses on the end surfaces of the beam was included in the formulation.

In this paper the shear stresses/strains that satisfied the shear stress free surface conditions were assumed to have their maximal values at the centroidal axis. This led to a displacement field that was fourth-order for beams with simple symmetrical cross-sections and third-order for beams with double symmetrical cross-sections. The equilibrium equations set on the basis of a vectorial approach led to the formulation of the governing equation (a fourth-order differential equation), efforts and deformations in terms of the transverse deflection. So closed-form solutions for efforts and deformations could be derived for various loading and support conditions, as well as element stiffness matrices.

2. Materials and Methods

2.1. First-Order Analysis of the Levinson Beam

The sign convention adopted for the loads, bending moments, shear forces, and displacements is illustrated in Figure 1

Let the displacements in x-, y-, and z- directions be denoted by u(x,y,z), v(x,y,z), and w(x,y,z), respectively. Observing the symmetry of the cross-section with respect to the z-axis, these displacements can be described as follows

$$u(x,y,z)=u(x,z),v(x,y,z)=0,w(x,y,z)=w(x)$$

(1)

The kinematic and constitutive relations are

$$\begin{array}{ll}{\epsilon }_x={\displaystyle \frac{\partial u}{\partial x}},& {\gamma }_{xz}={\displaystyle \frac{\partial u}{\partial z}}+{\displaystyle \frac{\partial w}{\partial x}}\\ {\sigma }_x=E{\epsilon }_x,& {\tau }_{xz}=G{\gamma }_{xz},\end{array}$$

(2)

respectively, where E and G are the Young’s modulus and the shear modulus, respectively.

In the Levinson Beam Theory (LBT), the displacement field is such that the shear stresses/strains vanish on the upper and lower surfaces of the beam. In addition, it is assumed in this study that the maximal shear stresses/strains appear in the centroidal axis.

Letting the height of the beam be h, and the z coordinates of the lower and upper surface be z_l = αh and z_u = - (1 - α) h, respectively, the displacement field can then be taken

$$\begin{array}{l}u(x,z)=u_0(x)+z\varphi \left(x\right)+\left(G\left(\alpha \right){\displaystyle \frac{z^4}{4h^3}}-F\left(\alpha \right){\displaystyle \frac{z^3}{3h^2}}\right)\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right),\begin{array}{cc}& w(x,z)=w(x),\end{array}\\ G\left(\alpha \right)={\displaystyle \frac{2\alpha -1}{{\alpha }^2{\left(1-\alpha \right)}^2}},\begin{array}{cc}& F\left(\alpha \right)={\displaystyle \frac{{\alpha }^3+{\left(1-\alpha \right)}^3}{{\alpha }^2{\left(1-\alpha \right)}^2}}\end{array},\end{array}$$

(3)

where u₀(x) is the displacement of the centroidal axis and φ(x) is the counterclockwise positive rotation of the cross-section at the centroidal axis. It is noted that for beams with double symmetrical cross sections α = 0.5 and the coefficients G(α) = 0, F(α) = 4, and so Equation (3) corresponds to the Levinson [2] solution for beams of rectangular cross sections. It follows the strains using Equation (2)

$$\begin{array}{l}{\epsilon }_x={\displaystyle \frac{d{u}_0(x)}{dx}}+z{\displaystyle \frac{d\varphi \left(x\right)}{dx}}+\left(G\left(\alpha \right){\displaystyle \frac{z^4}{4h^3}}-F\left(\alpha \right){\displaystyle \frac{z^3}{3h^2}}\right)\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right),\\ {\gamma }_{xz}=\left(1+G\left(\alpha \right){\displaystyle \frac{z^3}{h^3}}-F\left(\alpha \right){\displaystyle \frac{z^2}{h^2}}\right)\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right)\end{array}$$

(4)

The bending moment M(x), shear force V(x) and axial force N(x) are calculated using Equation (2) and (4) as follows

$$\begin{array}{l}M(x)={\displaystyle \underset{A}{\int }{\sigma }_{x}zdA}=EI{\displaystyle \frac{d\varphi \left(x\right)}{dx}}+E\left({\displaystyle \frac{G\left(\alpha \right)I_5}{4h^3}}-{\displaystyle \frac{F\left(\alpha \right)I_4}{3h^2}}\right)\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)\\ V(x)={\displaystyle \underset{A}{\int }{\tau }_{xz}dA}=G\left(A+{\displaystyle \frac{G\left(\alpha \right)I_3}{h^3}}-{\displaystyle \frac{F\left(\alpha \right)I}{h^2}}\right)\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right)\\ N(x)={\displaystyle \underset{A}{\int }{\sigma }_{x}dA}=EA{\displaystyle \frac{d{u}_0(x)}{dx}}+E\left({\displaystyle \frac{G\left(\alpha \right)I_4}{4h^3}}-{\displaystyle \frac{F\left(\alpha \right)I_3}{3h^2}}\right)\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)\end{array}$$

(5a–c)

where A and I are the area and the second moment of area of the cross-section, respectively, and I₃, I₄, and I₅ are defined as I_n = ∫zⁿdA with n = 3, 4, and 5. Let the following cross sectional values be

$$I_L=I+{\displaystyle \frac{G\left(\alpha \right)I_5}{4h^3}}-{\displaystyle \frac{F\left(\alpha \right)I_4}{3h^2}},\begin{array}{cc}& A_L=\end{array}A+{\displaystyle \frac{G\left(\alpha \right)I_3}{h^3}}-{\displaystyle \frac{F\left(\alpha \right)I}{h^2}},\begin{array}{cc}& S_L=\end{array}{\displaystyle \frac{G\left(\alpha \right)I_4}{4h^3}}-{\displaystyle \frac{F\left(\alpha \right)I_3}{3h^2}}$$

(6a-c)

Deriving both sides of Equation (5b) with respect to x and combining with (6b) yields

$${\displaystyle \frac{dV\left(x\right)}{dx}}=G{A}_L\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)$$

(7a)

Substituting Equation (7a) and (6a) into (5a) yields

$$M(x)=EI\left({\displaystyle \frac{1}{G{A}_L}}{\displaystyle \frac{dV\left(x\right)}{dx}}-{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)+E\left(I_L-I\right){\displaystyle \frac{1}{G{A}_L}}{\displaystyle \frac{dV\left(x\right)}{dx}}=-EI{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}+{\displaystyle \frac{E{I}_L}{G{A}_L}}{\displaystyle \frac{dV\left(x\right)}{dx}}$$

(7b)

Equation (7b) can be reformulated as follows to represent a moment − shear force − curvature relationship

$${\displaystyle \frac{d^{2}w(x)}{d{x}^2}}+{\displaystyle \frac{M(x)}{EI}}-{\displaystyle \frac{1}{G{A}_{L}I/I_L}}{\displaystyle \frac{dV(x)}{dx}}=0$$

(7c)

In the case of non-uniform heating, Equation (7c) becomes

$${\displaystyle \frac{d^{2}w(x)}{d{x}^2}}+{\displaystyle \frac{M(x)}{EI}}-{\displaystyle \frac{1}{G{A}_{L}I/I_L}}{\displaystyle \frac{dV(x)}{dx}}+{\displaystyle \frac{{\alpha }_T\Delta T}{h}}=0$$

(7d)

where α_T is the coefficient of thermal expansion, ΔT = ΔT_ls - ΔT_us is the difference between the temperature changes at lower surface (ΔT_ls) and upper surface (ΔT_us) of the beam, and h is as before the height of the beam.

The static equilibrium equations in first-order analysis are given by

$${\displaystyle \frac{dN(x)}{dx}}=0,{\displaystyle \frac{dV(x)}{dx}}=-p(x),V(x)={\displaystyle \frac{dM(x)}{dx}}$$

(8a-c)

Equation (8a-c) are developed as follows using Equation (5a-c) and (6a-c)

$$\begin{array}{l}EA{\displaystyle \frac{d^2{u}_0(x)}{d{x}^2}}+E{S}_L{\displaystyle \frac{d}{dx}}\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)=0\\ G{A}_L\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)=-p(x)\\ G{A}_L\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right)=EI{\displaystyle \frac{d^2\varphi \left(x\right)}{d{x}^2}}+E\left(I_L-I\right){\displaystyle \frac{d}{dx}}\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)\end{array}$$

(8d-f)

Let us formulate Equation (8e-f) for beams with rectangular cross-sections. Recalling that G(α) = 0 and F(α) = 4

$$\begin{array}{l}{A}_L=A-4I/h^2=bh-4{\displaystyle \frac{b{h}^3}{12}}{\displaystyle \frac{1}{h^2}}={\displaystyle \frac{2}{3}}A,I_4={\displaystyle \underset{-h/2}{\overset{h/2}{\int }}{z}^{4}dA}=b{\left.{\displaystyle \frac{z^5}{5}}\right|}_{-h/2}^{h/2}={\displaystyle \frac{b{h}^5}{80}}\\ I_L-I=-{\displaystyle \frac{F\left(\alpha \right)I_4}{3h^2}}=-4{\displaystyle \frac{b{h}^5}{80}}{\displaystyle \frac{1}{3h^2}}=-{\displaystyle \frac{1}{5}}I\end{array}$$

(8g)

Equation (8e-f) become

$$\begin{array}{l}{\displaystyle \frac{2}{3}}GA\left({\displaystyle \frac{d\varphi \left(x\right)}{dx}}+{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}\right)=-p(x)\\ {\displaystyle \frac{2}{3}}GA\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right)-{\displaystyle \frac{EI}{5}}\left(4{\displaystyle \frac{d^2\varphi \left(x\right)}{d{x}^2}}-{\displaystyle \frac{d^{3}w(x)}{d{x}^3}}\right)=0\end{array}$$

(8h-i)

Equation (8h-i) are the equilibrium equations derived by Levinson [2]. The bending moment can generally be expressed in terms of w(x) using Equation (7b) and (8b)

$$M(x)=-EI{\displaystyle \frac{d^{2}w(x)}{d{x}^2}}-{\displaystyle \frac{E{I}_L}{G{A}_L}}p(x)$$

(9)

The shear force can also be expressed in terms of w(x) using Equation (8c) and (9)

$$V(x)=-EI{\displaystyle \frac{d^{3}w(x)}{d{x}^3}}-{\displaystyle \frac{E{I}_L}{G{A}_L}}{\displaystyle \frac{dp(x)}{dx}}$$

(10)

Combining Equation (5b), (6b) and (10) yields

$$\begin{array}{l}G{A}_L\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right)=-EI{\displaystyle \frac{d^{3}w(x)}{d{x}^3}}-{\displaystyle \frac{E{I}_L}{G{A}_L}}{\displaystyle \frac{dp(x)}{dx}}\\ \to \varphi \left(x\right)=-{\displaystyle \frac{EI}{G{A}_L}}{\displaystyle \frac{d^{3}w(x)}{d{x}^3}}-{\displaystyle \frac{dw(x)}{dx}}-{\displaystyle \frac{E{I}_L}{{\left(G{A}_L\right)}^2}}{\displaystyle \frac{dp(x)}{dx}}\end{array}$$

(11)

Combining Equation (8b-c) and (9) yields the governing differential equation of the uniform Levinson beam with simply symmetrical cross section

$$EI{\displaystyle \frac{d^{4}w(x)}{d{x}^4}}=p(x)-{\displaystyle \frac{E{I}_L}{G{A}_L}}{\displaystyle \frac{d^{2}p(x)}{d{x}^2}}$$

(12)

The governing differential equation (Equation 12), the bending moment M (Equation 9), the shear force V (Equation 10), and the rotation of the cross section (Equation 11) are expressed in terms of w(x): the analytical solution can then be found.

For illustration we consider a beam of length l carrying a uniformly distributed load p. The transverse deflection w(x) can be expressed as follows using Equation (12)

$$EIw(x)={\displaystyle \frac{p}{24}}{x}^4+C_1{x}^3+C_2{x}^2+C_{3}x+C_4$$

(13)

It follows the bending moment, the shear force V, and the rotation of the cross section

$$\begin{array}{l}M(x)=-{\displaystyle \frac{p}{2}}{x}^2-6C_{1}x-2C_2-{\lambda }_L{\displaystyle \frac{I_L}{I}}p{l}^2,\\ V(x)=-px-6C_1,\\ EI\varphi \left(x\right)=-{\displaystyle \frac{p}{6}}{x}^3-3C_1{x}^2-\left(2C_2+{\lambda }_{L}p{l}^2\right)x-C_3-6{\lambda }_L{l}^2{C}_1,\\ {\lambda }_L={\displaystyle \frac{EI}{l^{2}G{A}_L}},\end{array}$$

(14a-d)

where λ_L is a bending shear factor. The integration constants C₁, C₂, C₃, and C₄ are then determined applying the boundary conditions expressed using Equation (13) and (14a-d). Results for some cases are presented in Appendix A.

It is noted that w(x) and φ(x) are determined independently from u₀(x). The latter is found from the following governing equation (8d) using (8e)

$$EA{\displaystyle \frac{d^2{u}_0(x)}{d{x}^2}}-{\displaystyle \frac{E{S}_L}{G{A}_L}}{\displaystyle \frac{dp\left(x\right)}{dx}}=0$$

(14e)

and the boundary conditions involving the axial force N(x) and u₀(x) at the beam’s ends.

Table 1 summarizes the Euler–Bernoulli beam equations, Timoshenko beam equations, and Levinson beam equations in terms of the transverse deflection w(x), whereby the rotation of cross-sections is counterclockwise positive and the Timoshenko shear coefficient is κ.

2.2. Stresses in the Cross-Section

In this section the stresses within the cross-section are expressed in terms of the axial force, the shear force, and the bending moment. They are calculated using Equation (2), (4), and (8e)

$$\begin{array}{l}{\sigma }_x=E{\displaystyle \frac{d{u}_0(x)}{dx}}+Ez{\displaystyle \frac{d\varphi \left(x\right)}{dx}}-{\displaystyle \frac{Ep\left(x\right)}{G{A}_L}}\left(G\left(\alpha \right){\displaystyle \frac{z^4}{4h^3}}-F\left(\alpha \right){\displaystyle \frac{z^3}{3h^2}}\right),\\ {\tau }_{xz}=G\left(1+G\left(\alpha \right){\displaystyle \frac{z^3}{h^3}}-F\left(\alpha \right){\displaystyle \frac{z^2}{h^2}}\right)\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right)\end{array}$$

(15a-b)

Substituting Equation (6c) and (8e) into (5c) yields

$$N(x)=EA{\displaystyle \frac{d{u}_0(x)}{dx}}-{\displaystyle \frac{E{S}_L}{G{A}_L}}p\left(x\right)$$

(16)

Equation (5a) can be reformulated as follows using (6a) and (8e) while (5b) is recalled using (6b)

$$\begin{array}{l}M(x)=EI{\displaystyle \frac{d\varphi \left(x\right)}{dx}}+E\left(I_L-I\right){\displaystyle \frac{-p\left(x\right)}{G{A}_L}}=EI{\displaystyle \frac{d\varphi \left(x\right)}{dx}}-{\displaystyle \frac{E\left(I_L-I\right)}{G{A}_L}}p\left(x\right)\\ V(x)=G{A}_L\left(\varphi \left(x\right)+{\displaystyle \frac{dw(x)}{dx}}\right)\end{array}$$

(17a-b)

From the equations above and reminding that E = 2G(1 + ν) it yields

$$\begin{array}{l}{\sigma }_x={\displaystyle \frac{N}{A}}+{\displaystyle \frac{M}{I}}z+2\left(1+\nu \right){\displaystyle \frac{p\left(x\right)}{A_L}}\left[{\displaystyle \frac{S_L}{A}}+\left({\displaystyle \frac{I_L}{I}}-1\right)z+F\left(\alpha \right){\displaystyle \frac{z^3}{3h^2}}-G\left(\alpha \right){\displaystyle \frac{z^4}{4h^3}}\right]\\ {\tau }_{xz}={\displaystyle \frac{V}{A_L}}\left(1+G\left(\alpha \right){\displaystyle \frac{z^3}{h^3}}-F\left(\alpha \right){\displaystyle \frac{z^2}{h^2}}\right)\end{array}$$

(18a-b)

In case of doubly symmetrical cross-sections G(α) = 0, F(α) = 4, S_L = 0, I_L = I – 4I₄/3h², A_L = A – 4I/h². Recalling that I₄ = ∫z⁴dA Equation (18a-b) become

$$\begin{array}{l}{\sigma }_x={\displaystyle \frac{N}{A}}+{\displaystyle \frac{M}{I}}z+{\displaystyle \frac{8\left(1+\nu \right)p\left(x\right)}{3A_L{h}^2}}z\left(z^2-{\displaystyle \frac{I_4}{I}}\right)\\ {\tau }_{xz}={\displaystyle \frac{V}{A_L}}\left(1-4{\displaystyle \frac{z^2}{h^2}}\right)\end{array}$$

(18c-d)

2.3. Element Stiffness Matrix

The sign convention for bending moments, shear forces, displacements, and rotation angles adopted for determining the element stiffness matrix in local coordinates is illustrated in Figure 2.

Let us define the element-end force vector and the element-end displacement vector as follows

$$\begin{array}{l}\overrightarrow{S}={\left[V_i{M}_i{V}_k{M}_k\right]}^T\\ \overrightarrow{V}={\left[w_i{\varphi }_i{w}_k{\varphi }_k\right]}^T\end{array}$$

(19a-b)

The element stiffness matrix in local coordinates of the Levinson beam is denoted by K_Ll. The relationship between the aforementioned vectors is as follows:

$$\overrightarrow{S}=K_{Ll}\times \overrightarrow{V}$$

(20)

Applying Equation (14a-d) with the distributed load p = 0 yields the following

$$\begin{array}{l}EIw(x)=C_1{x}^3+C_2{x}^2+C_{3}x+C_4\\ M(x)=-6C_{1}x-2C_2\\ V(x)=-6C_1\\ EI\varphi \left(x\right)=-3C_1{x}^2-2C_{2}x-C_3-6{\lambda }_L{l}^2{C}_1,\end{array}$$

(21a–d)

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 in combination with Equation (21a-d)

$$\begin{array}{l}{V}_i=-V(x=0)=6C_1\\ M_i=M(x=0)=-2C_2\\ V_k=V(x=l)=-6C_1\\ M_k=-M(x=l)=6C_{1}l+2C_2\end{array}$$

(22a–d)

Considering the sign conventions adopted for the displacements (see Figure 1) and rotations (counterclockwise positive) on the one hand and for displacements and rotations in the member stiffness matrix (see Figure 2) on the other hand, we can set following geometric compatibility boundary conditions in combination with Equation (21a-d)

$$\begin{array}{l}w(x=0)=w_i\to EI\times w_i=C_4\\ \varphi (x=0)=-{\varphi }_i\to EI\times {\varphi }_i=6{\lambda }_L{l}^2{C}_1+C_3\\ w(x=l)=w_k\to EI\times w_k=C_1{l}^3+C_2{l}^2+C_{3}l+C_4\\ \varphi (x=l)=-{\varphi }_k\to EI\times {\varphi }_k=3l^2\left(1+2{\lambda }_L\right)C_1+2C_{2}l+C_3\end{array}$$

(23a–d)

Combining Equation (20), (22a-d) and (23a-d) yields the 4 × 4 first-order element stiffness matrix of the Levinson beam

$$K_{Ll}=EI\times \left[\begin{array}{cccc}6& 0& 0& 0\\ 0& -2& 0& 0\\ -6& 0& 0& 0\\ 6l& 2& 0& 0\end{array}\right]\times {\left[\begin{array}{cccc}0& 0& 0& 1\\ 6{\lambda }_L{l}^2& 0& 1& 0\\ l^3& l^2& l& 1\\ 3l^2\left(1+2{\lambda }_L\right)& 2l& 1& 0\end{array}\right]}^{-1}$$

(24)

The development of Equation (24) yields

$$K_{Ll}={\displaystyle \frac{EI}{\left(1+{\Phi }_L\right)l^3}}\left[\begin{array}{cccc}12& 6l& -12& 6l\\ 6l& \left(4+{\Phi }_L\right)l^2& -6l& \left(2-{\Phi }_L\right)l^2\\ -12& -6l& 12& -6l\\ 6l& \left(2-{\Phi }_L\right)l^2& -6l& \left(4+{\Phi }_L\right)l^2\end{array}\right],\begin{array}{ccc}& & {\Phi }_L=12{\lambda }_L\end{array}$$

(25)

The element stiffness matrix (ESM) of the Levinson beam has the same formulation as that of the Timoshenko beam whereby the coefficients in the matrices are related as follows: A_LΦ_L = κAΦ_T. The ESM of the Euler−Bernoulli beam is obtained by setting Φ_L = 0.

Assuming the presence of a hinge at the right end, the sign convention for bending moments, shear forces, displacements, and rotations is illustrated in Figure 3.

The element-end force vector and the element-end displacement vector become

$$\begin{array}{l}\overrightarrow{S}={\left[V_i\begin{array}{cc}& M_i\begin{array}{cc}& V_k\end{array}\end{array}\right]}^T\\ \overrightarrow{V}={\left[w_i\begin{array}{cc}& {\varphi }_i\begin{array}{cc}& w_k\end{array}\end{array}\right]}^T\end{array}$$

(26)

The analysis is done similarly and the 3 × 3 element stiffness matrix of the Levinson beam is then

$$K_{Ll}={\displaystyle \frac{EI}{\left(1+3{\lambda }_L\right)l^3}}\left[\begin{array}{ccc}3& 3l& -3\\ 3l& 3l^2& -3l\\ -3& -3l& 3\end{array}\right]$$

(27)

As before the element stiffness matrix of the Levinson beam has the same formulation as that of the Timoshenko beam whereby the coefficients in the matrices are related as follows: A_Lλ_L = κAλ_T.

2.4. Strong and Weak Formulation of the Boundary Conditions

The analytical solution presented in this paper satisfies the boundary conditions (BC’s) through the stress resultants (axial force, shear force and bending moment) and the rotation of cross-section only at the centroidal axis: this represents a weak formulation of the BC’s. Strong formulations of the BC’s expressed in terms of stresses and displacements are listed below for different edge conditions whereby the displacement component u₀(x) is ignored since it does not contribute to M(x) and V(x).

Clamped Edge

From the strong formulation the displacements u(x, z) and w(x) are zero at any position z. From Equation (3) it yields φ = 0, dw/dx = 0, w = 0. Consequently the shear stresses and the shear force vanish regarding Equation (4) and (5b), this representing an inconsistency. By applying a weak formulation of the BC’s (φ = 0, w = 0) this inconsistency disappears.

Free Edge

The strong formulation implies the normal and shear stresses σ and τ are zero at any position z. Equation (15a-b) yields dφ/dx = 0, p(x) = 0, φ + dw/dx = 0. Hence, the absence of distributed loading is a condition for the satisfaction in a strong sense of the BC’s. The weak formulation is M = 0, V = 0.

Simply Supported Edge

From the strong formulation the normal stresses σ = 0 and the transverse deflection w = 0 at any position z. Equation (15a-b) implies dφ/dx = 0, p(x) = 0, w = 0. The absence of distributed loading is a condition for the satisfaction in a strong sense of the BC’s. The weak formulation is M = 0, w = 0.

In summary, on the one hand it is observed that the strong formulation of the BC’s is characterized with three conditions: this required a beam model described with a sixth-order governing differential equation. On the other hand the weak formulation is characterized with two conditions at each edge, which is consistent with the fourth-order governing differential equation of the Levinson beam model.

2.5. Inconsistencies in the Determination of Shear Stresses

In this study the shear stresses were determined using the stress-strain relationship or constitutive law (Equation (2d)). However they can be calculated using the equilibrium-method as it is done in the Euler−Bernoulli beam theory. Here a beam element obtained by performing plane cuts at x, x + dx, and z is considered. The equilibrium equation in x direction is set on this element using Equation (15a-b) for stresses, i.e., the shear stresses τ_zx = τ_xz have to balance-out the variation of the axial force apart from z as follows

$$\begin{array}{l}b\left(z\right){\tau }_{xz}={\displaystyle \underset{z}{\overset{z_l}{\int }}\frac{\partial {\sigma }_x}{\partial x}dA}\\ \begin{array}{ccc}& & \end{array}={\displaystyle \underset{z}{\overset{z_l}{\int }}\frac{\partial }{\partial x}\left[E\frac{d{u}_0(x)}{dx}+Ez\frac{d\varphi \left(x\right)}{dx}-\frac{Ep\left(x\right)}{G{A}_L}\left(G\left(\alpha \right)\frac{z^4}{4h^3}-F\left(\alpha \right)\frac{z^3}{3h^2}\right)\right]dzdy}\end{array}$$

(28)

The term on the right-hand side of Equation (28) is a fifth-order polynomial with respect to z while τ_xz (Equation (15b)) is a third-order polynomial with respect to z.

Therefore, similarly to Euler−Bernoulli beam theory and Timoshenko beam theory, an inconsistency remains between the shear stresses that are obtained based on the equilibrium-method and the shear stresses obtained by the constitutive law.

3. Results and Discussion

3.1. Beams with Rectangular Cross-Section

Normal and Shear Stresses

Equation (18c-d) presented the stresses for doubly symmetrical cross-sections. For rectangular cross-sections A_L = 2/3A and I₄ = bh⁵/80 were calculated in Equation (8g): it follows the stresses

$$\begin{array}{l}{\sigma }_x={\displaystyle \frac{N}{A}}+{\displaystyle \frac{M}{I}}z+{\displaystyle \frac{8\left(1+\nu \right)p\left(x\right)}{3A_L{h}^2}}z\left(z^2-{\displaystyle \frac{I_4}{I}}\right)={\displaystyle \frac{N}{A}}+{\displaystyle \frac{M}{I}}z+{\displaystyle \frac{\left(1+\nu \right)p}{60I}}z\left(20z^2-3h^2\right)\\ {\tau }_{xz}={\displaystyle \frac{V}{A_L}}\left(1-4{\displaystyle \frac{z^2}{h^2}}\right)={\displaystyle \frac{V}{8I}}\left(h^2-4z^2\right)\end{array}$$

(29)

These results correspond to those obtained by Karttunen [9]

Element Stiffness Matrix

Equation (25) presented the 4 × 4 element stiffness matrix with the parameter Φ_L = 12λ_L. For rectangular cross-section

$${\Phi }_L=12{\lambda }_L=12{\displaystyle \frac{EI}{l^{2}G{A}_L}}=12\times 2\left(1+\nu \right){\displaystyle \frac{b{h}^3}{12}}{\displaystyle \frac{1}{l^2}}{\displaystyle \frac{1}{\frac{2}{3}bh}}=3\left(1+\nu \right){\displaystyle \frac{h^2}{l^2}}$$

(30)

This result corresponds to that obtained by Karttunen [9]

Deflection Curve of Beam of Length 2l Simply Supported at Both Ends

A beam of uniform cross-section and length 2l simply supported at both ends and carrying a uniformly distributed load p is considered. The deflection curve is to determine whereby the origin is taken at the middle for comparison purpose with Levinson [2]. The boundary conditions are M(±l) = 0 and w(±l) = 0 at both ends while the expressions of M(x) and w(x) were given in Equation (13) and (14a).

The solution for the deflection curve of a beam of length 2l with an arbitrary simply symmetrical cross-section is

$$EIw(x)={\displaystyle \frac{p{l}^4}{24}}\left({\displaystyle \frac{x^4}{l^4}}-6{\displaystyle \frac{x^2}{l^2}}+5\right)+2{\lambda }_L{\displaystyle \frac{I_L}{I}}p{l}^4\left(1-{\displaystyle \frac{x^2}{l^2}}\right)$$

(31)

For a beam of length 2l having a rectangular cross-section Equation (8g) yields I_L = 4/5I and

$$\begin{array}{l}{\lambda }_L{\displaystyle \frac{I_L}{I}}={\displaystyle \frac{EI}{{\left(2l\right)}^{2}G{A}_L}}{\displaystyle \frac{4}{5}}=2\left(1+\nu \right){\displaystyle \frac{b{h}^3}{12}}{\displaystyle \frac{1}{4l^2}}{\displaystyle \frac{1}{\frac{2}{3}bh}}{\displaystyle \frac{4}{5}}=\left(1+\nu \right){\displaystyle \frac{h^2}{20l^2}}\\ \to EIw(x)={\displaystyle \frac{p{l}^4}{24}}\left({\displaystyle \frac{x^4}{l^4}}-6{\displaystyle \frac{x^2}{l^2}}+5\right)+{\displaystyle \frac{1}{10}}p{l}^4\left(1+\nu \right){\displaystyle \frac{h^2}{l^2}}\left(1-{\displaystyle \frac{x^2}{l^2}}\right)\end{array}$$

(32)

This result corresponds to that obtained by Levinson [2]

4. Conclusions

The Levinson beam theory (LBT) is treated by many authors for beams with rectangular cross-sections. In this paper LBT was extended to simple symmetrical beams. A displacement field which includes warping of the cross-section and satisfies the shear free conditions on the lower and upper surfaces was presented and two coupled differential equations combining the transverse deflection and the rotation of the cross-section at the centroidal axis were derived. The analysis was simplified by the formulation of the governing equation, the efforts and deformations in terms of the transverse deflection and so closed-form solutions were presented for various loading and support conditions, as well as element stiffness matrices. It was shown that a strong formulation of the boundary conditions (stress and displacement related) presented many inconsistencies while the weak formulation (related to the efforts and the rotation of the cross-section at the centroidal axis) utilized in this study gave results that were in agreement with those in the literature for beams with rectangular cross-sections. However, an inconsistency remained between the shear stresses that are obtained based on the equilibrium-method and the shear stresses obtained by the constitutive law.

The following aspects of LBT not addressed in this study could be examined in future research:

• ✓ Second-order analysis of beams
• ✓ Buckling analysis of beams
• ✓ Vibration analysis of beams