Isotropic elastic beam structures

7.3. Isotropic elastic beam structures#

In this section, theorems for elastic structures are specialized for beam structures.

todo Beam model used here…

7.3.1. Strong formulation of the problem#

Indefinite equilibrium. External distributed loads are equilibrated by internal actions, resultant of stress field on beam sections.

\[\begin{split}\begin{aligned} 0 & = F_z'(z) + n(z) && \text{(axial loads)} \\ 0 & = F_x'(z) + f_x(z) && \text{(shear loads)} \\ 0 & = F_y'(z) + f_y(z) && \\ 0 & = M_x'(z) - F_y(z) + m_x(z) && \text{(bending)} \\ 0 & = M_y'(z) + F_x(z) + m_y(z) && \\ 0 & = M_z'(z) + m_z(z) && \text{(torsion)} \\ \end{aligned}\end{split}\]

Kinematics.

Constitutive equations. Under the assumptions … (kinematic assumptions, decoupling,…),

\[\begin{split}\begin{aligned} F_z(z) & = EA s'_z(z) \\ F_x(z) & = \chi^{-1}_x GA ( s'_x(z) - \theta_y(z) ) \\ F_y(z) & = \chi^{-1}_y GA ( s'_y(z) + \theta_x(z) ) \\ M_x(z) & = EJ_x \theta'_x(z) \\ M_y(z) & = EJ_y \theta'_y(z) \\ M_z(z) & = GJ_y \theta'_z(z) \\ \end{aligned}\end{split}\]

With thermal strains due to temperature distribution $T(z) = T_0(z) + \Delta_x T(z) \frac{x}{h_x} + \Delta_y T(z) \frac{y}{h_y}$, these constitutive equations hold for the mechanical part only, while the most general constitutive equations read

\[\begin{split}\begin{aligned} u'_z & = u^{mech \ '}_z + u^{th \ '}_z = \dfrac{T_z}{EA} + \alpha \Delta T \\ \theta'_x & = \theta^{mech \ '}_x + \theta^{th \ '}_x = \dfrac{M_x}{EJ_x} + \alpha \dfrac{\Delta_y T}{h_y} \\ \theta'_y & = \theta^{mech \ '}_y + \theta^{th \ '}_y = \dfrac{M_y}{EJ_y} - \alpha \dfrac{\Delta_x T}{h_x} \\ ... \\ \end{aligned}\end{split}\]

Bernoulli beam. If Bernoulli kinematic assumption

\[\begin{split}\begin{aligned} \theta_x(z) & = - u'_y(z) \\ \theta_y(z) & = u'_x(z) \\ \end{aligned}\end{split}\]

holds, bending equilibrium equations read

\[\begin{split}\begin{aligned} 0 & = f_x + T_x' = f_x - M''_y = f_x - ( EJ_y \theta'_y )'' = f_x - ( EJ_y u''_x )'' \\ 0 & = f_y + T_y' = f_y + M''_x = f_y + ( EJ_x \theta'_x )'' = f_y - ( EJ_x u''_y )'' \\ \end{aligned}\end{split}\]

7.3.2. Weak formulations of the problem#

7.3.2.1. Weak formulation of equilibrium conditions#

Under the assumption of negligible contribution of shear stress and deformation, and Bernoulli kinematic assumption, the weak form of equilibrium equation is derived as follows: axial and bending indefinite equilibrium equations for each beam (or structural element in general) are multiplied by arbitrary test functions, these products are integrated over the beam length. Then, natural boundary conditions are applied.

Example 7.1 (Finite element method. Example: 1-element clamped Bernoulli beam)

Approximation 1. Choosing only one test function $u_{y,2}(z) = z^2$ (that satisfies the essential boundary conditions, ($u_{y,2}(0)=0$, $u'_{y,2}(0)=0$), and assuming that the $y$-displacement $s_y$ is proportional to that test function (same test and base function, to get a symmetric problem - even when multi-dimensional),

\[s_y(z) = a_2 z^2\]

and the bending moment reads

\[M_x(z) = - EJ_x s_y''(z) = - 2 EJ_x a_2 \ ,\]

the weak formulation of equilibrium equation reads (assuming Bernoulli beam, no distributed loads, and lumped force $F_y$ at the free end as the only load)

\[\begin{split}\begin{aligned} 0 & = - \int_{z=0}^{b} ( 2 EJ_x a_2 ) ( 2 ) \, dz + F_y b^2 = \\ & = - 4 EJ_x b \, a_2 + F_y b^2 \ , \end{aligned}\end{split}\]

so that the only coefficient of the approximation is $a_2 = \frac{b}{4 EJ} F_y$. Thus the approximate solution of the problem is $s_y(z) = \frac{F_y b}{4 EJ_x} z^2$, and the displacement of the free end is $w(b) = \frac{Fb^3}{4 EJ_x}$.

Comment. The displacement of the approximate solution is smaller than the displacement of the exact solution $w_B = \frac{Fb^3}{3 EJ_x}$, as the approximated solution can’t represent the exact solution, $s_y(z) = \frac{Fb^3}{EJ}\left( - \frac{1}{6}\left(\frac{z}{b}\right)^3 + \frac{1}{2}\left(\frac{z}{b}\right)^2 \right)$, since it doesn’t contain any $z^3$ term, and this introduce an extra numerical additional constraint that makes the approximate model stiffer than the continuous model.

Approximation 2.

\[\begin{split}s_y(z) = a_2 z^2 + a_3 z^3 = \begin{bmatrix} z^2 & z^3 \end{bmatrix} \begin{bmatrix} a_2 \\ a_3 \end{bmatrix}\end{split}\]

\[\begin{split}s''_y(z) = 2 a_2 + 6 a_3 z = \begin{bmatrix} 2 & 6 z \end{bmatrix} \begin{bmatrix} a_2 \\ a_3 \end{bmatrix}\end{split}\]

\[\begin{split}\begin{aligned} \mathbf{0} & = - \int_{z=0}^{b} \begin{bmatrix} 2 \\ 6 z \end{bmatrix} EJ_x \begin{bmatrix} 2 & 6 z \end{bmatrix} \, dz \, \begin{bmatrix} a_2 \\ a_3 \end{bmatrix} + \begin{bmatrix} b^2 \\ b^3 \end{bmatrix} F_y = \\ & = - EJ_x \begin{bmatrix} 4 b & 6 b^2 \\ 6 b^2 & 12 b^3 \end{bmatrix} \, \begin{bmatrix} a_2 \\ a_3 \end{bmatrix} + \begin{bmatrix} b^2 \\ b^3 \end{bmatrix} F_y \ , \end{aligned}\end{split}\]

and the solution of the system provides the value of the coefficients of the linear combination in the approximate solution

\[\begin{split}\begin{aligned} \begin{bmatrix} a_2 \\ a_3 \end{bmatrix} = \dfrac{1}{12 b^4} \begin{bmatrix} 12 b^3 & - 6 b^2 \\ - 6 b^2 & 4 b \end{bmatrix} \begin{bmatrix} b^2 \\ b^3 \end{bmatrix} \dfrac{F_y}{E J_x} = \begin{bmatrix} \frac{b}{2} \\ -\frac{1}{6} \end{bmatrix} \dfrac{F_y}{EJ_x} \ , \end{aligned}\end{split}\]

This approximation allows to get the exact solution of the continuous problem,

\[s_y(z) = \frac{F_y b}{2 EJ_x} z^2 - \frac{F_y}{6 EJ_x} z^3 \ .\]

7.3.2.2. Weak formulation of congruence conditions#

Timoshenko beam

Neglecting transverse strain and warping contributions (their contribution is zero when multiplied by constant and linear functions and integrated), strain field in a Timoshenko beam reads

\[\begin{split}\begin{aligned} \varepsilon_{zz} & = s'_{Pz} + y \theta'_x - x \theta'_y \\ \varepsilon_{xx} & = 0 \\ \varepsilon_{yy} & = 0 \\ 2 \varepsilon_{zx} & = s'_{Px} - \theta_y - y \theta'_z \\ 2 \varepsilon_{zy} & = s'_{Py} + \theta_x + x \theta'_z \\ 2 \varepsilon_{xy} & = 0 \\ \end{aligned}\end{split}\]

A weak form of congruence conditions is obtained multiplying by test functions…(todo discuss the choice of test functions)

\[\begin{split}\begin{aligned} 0 & = - \int_V \left\{ (\Sigma_z + y \Phi_{x} - x \Phi_{y}) (- \varepsilon_{zz} + s'_{Pz} + y \theta'_x - x \theta'_y ) + ( \Sigma_{x} - y \Phi_z ) ( - 2\varepsilon_{zx} + s'_{Px} - \theta_y - y \theta'_z ) + ( \Sigma_{y} - x \Phi_z ) ( - 2\varepsilon_{zy} + s'_{Py} + \theta_x - x \theta'_z ) \right \} = \\ & = \int_V \left\{ ( \Sigma_z + y \Phi_{x} - x \Phi_{y} ) \varepsilon_{zz} + ( 2 \Sigma_x ) \varepsilon_{zx} + ( 2 \Sigma_y ) \varepsilon_{zy} \right\} + \\ & \qquad + \int_V \left\{ ( \Sigma_z + y \Phi_{x} - x \Phi_{y} )' ( s_{Pz} + y \theta_x - x \theta_y ) + ( \Sigma_x' - y \Phi'_z ) ( s_{Px} - y \theta_z ) + ( \Sigma_x - y \Phi_z ) \theta_y + ( \Sigma_y' - x \Phi'_z ) ( s_{Py} - x \theta_z ) - ( \Sigma_y - x \Phi_z ) \theta_x \right\} + \\ & \qquad - \left[ \int_A \left\{ ( \Sigma_z + y \Phi_{x} - x \Phi_{y} ) ( s_{Pz} + y \theta_x - x \theta_y ) + ( \Sigma_x - y \Phi_z ) (s_{Px} - y \theta_z) + ( \Sigma_y - z \Phi_z ) ( s_{Py} - x \theta_z ) \right\} \right]_{\partial l} \end{aligned}\end{split}\]

For simplicity, considering first a structurally decoupled system, so that static moments are zero and the inertia tensor is diagonal,

\[\begin{split}\begin{aligned} 0 & = \int_V \left\{ ( \Sigma_z + y \Phi_{x} - x \Phi_{y} ) \varepsilon_{zz} + 2 ( \Sigma_x - y \Phi_z ) \varepsilon_{zx} + 2 ( \Sigma_y - x \Phi_z ) \varepsilon_{zy} \right\} + \\ & \qquad + \int_\ell \left\{ \Sigma'_x A s_{Px} + \Sigma'_y A s_{Py} + \Sigma'_z A s_{Pz} + \left( \Phi'_x J_x - \Sigma_y A \right) \theta_x + ( \Phi'_y J_y + \Sigma_x A ) \theta_y + \Phi'_z J_z \theta_z \right\} + \\ & \qquad - \left[ \Sigma_z A s_{Pz} + \Phi_x J_x \theta_x + \Phi_y J_y \theta_y + \Sigma_{x} s_{Px} + \Sigma_y s_{Py} + \Phi_z \theta_z \right]_{\partial l} \end{aligned}\end{split}\]

If the test functions are related to equilibrated internal actions,

\[\begin{split}\begin{aligned} 0 & = \widetilde{F}_z'(z) + \widetilde{f}_z(z) && \text{(axial loads)} \\ 0 & = \widetilde{F}_x'(z) + \widetilde{f}_x(z) && \text{(shear loads)} \\ 0 & = \widetilde{F}_y'(z) + \widetilde{f}_y(z) && \\ 0 & = \widetilde{M}_x'(z) - \widetilde{F}_y(z) + \widetilde{m}_x(z) && \text{(bending)} \\ 0 & = \widetilde{M}_y'(z) + \widetilde{F}_x(z) + \widetilde{m}_y(z) && \\ 0 & = \widetilde{M}_z'(z) + \widetilde{m}_z(z) && \text{(torsion)} \\ \end{aligned}\end{split}\]

and defined as $\widetilde{F}_i = \Sigma_i A$, $\widetilde{M}_i = \Phi_i J_i$, the weak form of the problem becomes

\[\begin{split}\begin{aligned} 0 & = \int_V \left\{ \left( \frac{\widetilde{F}_z}{A} + y \frac{\widetilde{M}_x}{J_x} - x \frac{\widetilde{M}_y}{J_y} \right) \varepsilon_{zz} + 2 \left( \dfrac{\widetilde{F}_x}{A} - y \dfrac{\widetilde{M}_z}{J_z} \right) \varepsilon_{zx} + 2 \left( \dfrac{\widetilde{F}_y}{A} - x \dfrac{\widetilde{M}_z}{J_z} \right) \varepsilon_{zy} \right\} + \\ & \qquad - \int_\ell \left\{ \widetilde{f}_x s_{Px} + \widetilde{f}_y s_{Py} + \widetilde{f}_z s_{Pz} + \widetilde{m}_x \theta_x - \widetilde{m}_y \theta_y - \widetilde{m}_z \theta_z \right\} + \\ & \qquad - \left[ \widetilde{F}_x s_{Px} + \widetilde{F}_y s_{Py} + \widetilde{F}_z s_{Pz} + \widetilde{M}_x \theta_x + \widetilde{M}_y \theta_y + \widetilde{M}_z \theta_z \right]_{\partial l} \ . \end{aligned}\end{split}\]

7.3.2.3. Principle of virtual work#

Starting from the weak form of equilibrium conditions, and choosing the test functions $u_i = \delta s_i$, $\phi_i = \delta \theta_i$ to be variations of congruent displacement fields,

\[\begin{split}\begin{aligned} & \delta s_x = 0 && \text{where $x$-transverse displacement is prescribed} \\ & \delta s_y = 0 && \text{where $y$-transverse displacement is prescribed} \\ & \delta s_z = 0 && \text{where $z$-axial displacement is prescribed} \\ & \delta \theta_x = 0 && \text{where $x$-rotation is prescribed} \\ & \delta \theta_y = 0 && \text{where $y$-rotation is prescribed} \\ & \delta \theta_z = 0 && \text{where $z$-rotation is prescribed} \\ \end{aligned}\end{split}\]

7.3.2.4. Principle of complementary virtual work#

Timoshenko beam

If the test functions are related to equilibrated internal actions,

\[\begin{split}\begin{aligned} 0 & = \delta \widetilde{F}_z'(z) && \text{(axial loads)} \\ 0 & = \delta \widetilde{F}_x'(z) && \text{(shear loads)} \\ 0 & = \delta \widetilde{F}_y'(z) && \\ 0 & = \delta \widetilde{M}_x'(z) - \delta \widetilde{F}_y(z) && \text{(bending)} \\ 0 & = \delta \widetilde{M}_y'(z) + \delta \widetilde{F}_x(z) && \\ 0 & = \delta \widetilde{M}_z'(z) && \text{(torsion)} \\ \end{aligned}\end{split}\]

and defined as $\widetilde{F}_i = \Sigma_i A$, $\widetilde{M}_i = \Phi_i J_i$, the weak form of the problem becomes

\[\begin{split}\begin{aligned} 0 & = \int_V \left\{ \left( \frac{\delta \widetilde{F}_z}{A} + y \frac{\delta \widetilde{M}_x}{J_x} - x \frac{\delta \widetilde{M}_y}{J_y} \right) \varepsilon_{zz} + 2 \left( \dfrac{\delta \widetilde{F}_x}{A} - y \dfrac{\delta \widetilde{M}_z}{J_z} \right) \varepsilon_{zx} + 2 \left( \dfrac{\delta \widetilde{F}_y}{A} - x \dfrac{\delta \widetilde{M}_z}{J_z} \right) \varepsilon_{zy} \right\} + \\ & \qquad - \left[ \delta \widetilde{F}_x s_{Px} + \delta \widetilde{F}_y s_{Py} + \delta \widetilde{F}_z s_{Pz} + \delta \widetilde{M}_x \theta_x + \delta \widetilde{M}_y \theta_y + \delta \widetilde{M}_z \theta_z \right]_{\partial l / S_D} \ . \end{aligned}\end{split}\]

Introducing the constitutive law for an isotropic elastic beam with structural decoupling, it’s possible to write strain as a function of the internal actions

\[\begin{split}\begin{aligned} \varepsilon_{zz} & = \dfrac{F_z}{EA} + y \dfrac{M_x}{J_x} - x \dfrac{M_y}{J_y} \\ 2 \varepsilon_{zx} & = \chi_x \dfrac{F_x}{GA} \\ 2 \varepsilon_{zy} & = \chi_y \dfrac{F_y}{GA} \ , \end{aligned}\end{split}\]

s.t. the PCVW (remember structural decoupling) becomes

\[\begin{split}\begin{aligned} 0 & = \int_{\ell} \left\{ \frac{\delta \widetilde{F}_z F_z}{E A} + \frac{\delta \widetilde{F}_x F_x}{\chi^{-1}_x G A} + \frac{\delta \widetilde{F}_y F_y}{\chi^{-1}_y G A} + \frac{\delta \widetilde{M}_x M_x}{EJ_x} + \frac{\delta \widetilde{M}_y M_y}{EJ_y} + \frac{\delta \widetilde{M}_z M_z}{GJ_z} \right\} + \\ & \qquad - \left[ \delta \widetilde{F}_x s_{Px} + \delta \widetilde{F}_y s_{Py} + \delta \widetilde{F}_z s_{Pz} + \delta \widetilde{M}_x \theta_x + \delta \widetilde{M}_y \theta_y + \delta \widetilde{M}_z \theta_z \right]_{\partial l / S_D} \ . \end{aligned}\end{split}\]

7.3.2.5. Principle of stationariety of total potential energy#

7.3.2.6. Principle of stationariety of total complementary potential energy#

Example 7.2 (Clamp-cart beam)

A beam of length $b$ is clamped in $A$ and constrained with a cart in $B$ preventing transverse displacement. A transverse distributed uniform load $q$ is applied. Thermal deformation is induced by a linear distribution of $\Delta T$ across all the sections of the beam. Assuming slender beam model and considering only the bending deformation, solve the (hyperstatic) structure: determine the internal bending moment $M(z)$ and the displacement $w(z)$.

Let’s solve the problem with different approaches, after calling $X$ the cart reaction and choosing it as the independent unknown hyperstatics. The reaction in $A$ are $V_A = - X + qb$, $H_A = 0$, $M_A = - \frac{qb^2}{2} + Xb$ (clockwise), and the internal bending moment from $A: \, z=0$ to $B: \, z=b$ (counter-clockwise)

\[\begin{split}\begin{aligned} M(z) & = M_A + V_A \, z - \frac{q z^2}{2} = \\ & = -\frac{qb^2}{2} + Xb + \left( - X + qb \right) z - \frac{q z^2}{2} = \\ & = -\frac{1}{2} q \left( b - z \right)^2 + X ( b - z ) \ . \end{aligned}\end{split}\]

Method 1 - Elastic line. The axial equilibrium is trivial with axial internal action $N(z) = 0$, and displacement $u(z) = 0$. Focusing on bending equilibrium and transverse displacement

\[\theta'(z) = \theta^{mech \ '}(z) + \theta^{th \ '}(z) = \dfrac{M(z)}{EJ} + \dfrac{\alpha \Delta T}{h} \ ,\]

with the approximation $\theta(z) \simeq w'(z)$, the problem is governed by the differential problem

\[ w''(z) = \frac{1}{EJ} \left[ -\frac{1}{2} q \left( b - z \right)^2 + X ( b - z ) \right] + \frac{\alpha \Delta T}{h} \ , \]

supplied with the boundary conditions $w(z=0) = 0$, $w'(z=0) = 0$, $w(z=b) = 0$. The solution reads

\[\begin{split}\begin{aligned} w'(z) & = \frac{q}{EJ} \left[ - \frac{b^2 z}{2} + \frac{b z^2}{2} - \frac{z^3}{6} \right] + \frac{X}{EJ} \left[ \left( bz - \frac{z^2}{2} \right) \right] + \frac{\alpha \Delta T}{h}z + A \\ w (z) & = \frac{q}{EJ} \left[ - \frac{b^2 z^2}{4} + \frac{b z^3}{6} - \frac{z^4}{24} \right] + \frac{X}{EJ} \left[ \left( b \frac{z^2}{2} - \frac{z^3}{6} \right) \right] + \frac{\alpha \Delta T}{2h}z^2 + Az + B \ . \end{aligned}\end{split}\]

Boundary conditions in $z=0$ forces $A = 0$, $B = 0$. The boundary condition in $B: \, z = b$ is used to evaluate the hyperstatics $X$

\[0 = w(b) = -\frac{qb^4}{8 EJ} + \frac{Xb^3}{3EJ} + \frac{\alpha \Delta T \, b^2}{2 h} \ ,\]

so that:

the hyperstatics is

\[X = \frac{3}{8} qb - \frac{3}{2}\frac{\alpha \Delta T \, EJ}{hb} \ ,\]
the internal bending moment is

\[M(z) = - \frac{qz^2}{2} + \left( \frac{5}{8} qb + \frac{3}{2} \frac{\alpha \Delta T \, EJ}{h} \right) z + \left( - \frac{1}{8} q b^2 - \frac{3}{2} \alpha \Delta T \, EJ b \right) \]
the transverse displacement is

\[w(z) = -\frac{q z^4}{24 \, EJ} + \left( \frac{5}{8} qb + \frac{3}{2} \frac{\alpha \Delta T \, EJ}{h} \right) \frac{z^3}{6} + \left( - \frac{1}{8} q b^2 - \frac{3}{2} \alpha \Delta T \, EJ b \right) \frac{z^2}{2} \ .\]

Proof, $\ w(b) = 0$. $-\frac{1}{24} + \frac{5}{48} - \frac{1}{16} = \frac{-2+5-3}{48} = 0 .$

Method 2 - Force method. An equilibrated solution with unit external loads at the hyperstatics, $\widetilde{X} = 1$, has intenral bending moment $\widetilde{M}(z) = X ( b - z)$. The PCVW with this equilibrated solution as the test function in the weak formulation reads

\[\begin{split}\begin{aligned} 0 & = \int_{z=0}^{b} \widetilde{M}(z) \theta'(z) \, dz - \widetilde{X} \underbrace{w_B}_{=0} = \\ & = \int_{z=0}^{b} \widetilde{M}(z) \left( \frac{M(z)}{EJ} + \frac{\alpha \Delta T}{h} \right) \, dz = \\ & = \int_{z=0}^{b} \left( b - z \right) \left[ \frac{1}{EJ} \left( -\frac{1}{2} q (b-z)^2 + X(b-z) \right) + \frac{\alpha \Delta T}{h} \right] \, dz = \\ & = - \dfrac{1}{8} \frac{qb^4}{EJ} + \frac{1}{3}\frac{Xb^3}{EJ} + \frac{\alpha \Delta T \, b^2}{2 h} \ , \end{aligned}\end{split}\]

and thus the same expression of the hyperstatics is found, $X = \frac{3}{8} qb - \frac{3}{2} \frac{\alpha \Delta T \, EJ}{hb}$.

Method 3 - Stationariety of total complementary potential energy. Total complementary potential energy reads

\[\Pi^* = \int_{z=0}^{b} \frac{1}{2} \frac{M^2(z; X)}{EJ} \, dz + \int_{z=0}^{b} M(z; X) \frac{\alpha \Delta T}{h} \, dz \ ,\]

and its derivative1 w.r.t. $X$ (the only independent variable; equilibrium and constitutive law are already used to get the latter expression),

\[\begin{aligned} 0 & = \frac{\partial \Pi^*}{\partial X} = \int_{z=0}^{b} \frac{M(z; X)}{EJ} \frac{\partial M}{\partial X}(z; X) \, dz + \int_{z=0}^{b} \frac{\partial M}{\partial X}(z;X) \frac{\alpha \Delta T}{h} \, dz \ . \end{aligned}\]

Using the expression of teh internal bending moment $M(z; X) = -\frac{1}{2} q (b-z)^2 + X(b-z)$, the same expression as the one provided by method 2 - force method is found. Thus, doing algebra properly, this method gives the same results as the other two methods.

1: The variation of bending moment becomes $\delta M(z; X) = \frac{\partial M}{\partial X} \delta X$, and thus $\delta \Pi^* = \frac{\partial \Pi^*}{\partial X} \delta X$. If the total complementary potential energy functional needs to be stationary, $\delta \Pi^*$ for every possible $\delta X$, it follows $\frac{\partial \Pi^*}{\partial X} = 0$.