# Nonlinear finite elements/Euler Bernoulli beams

## Euler-Bernoulli Beam

### Displacements

{\begin{aligned}u_{1}&=u_{0}(x)-z{\cfrac {dw_{0}}{dx}}\\u_{2}&=0\\u_{3}&=w_{0}(x)\end{aligned}} ### Strains

$\varepsilon _{11}=\varepsilon _{xx}=\varepsilon _{xx}^{0}+z\varepsilon _{xx}^{1}$ {\begin{aligned}\varepsilon _{xx}^{0}&={\cfrac {du_{0}}{dx}}+{\frac {1}{2}}\left({\cfrac {dw_{0}}{dx}}\right)^{2}\\\varepsilon _{xx}^{1}&=-{\cfrac {d^{2}w_{0}}{dx^{2}}}\end{aligned}} ## Strain-Displacement Relations

$\varepsilon _{ij}={\frac {1}{2}}\left({\frac {\partial u_{i}}{\partial x_{j}}}+{\frac {\partial u_{j}}{\partial x_{i}}}\right)+{\frac {1}{2}}\left({\frac {\partial u_{m}}{\partial x_{i}}}{\frac {\partial u_{m}}{\partial x_{j}}}\right)$ {\begin{aligned}\varepsilon _{11}&={\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{1}}}+{\frac {\partial u_{1}}{\partial x_{1}}}\right)+{\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{1}}}{\frac {\partial u_{1}}{\partial x_{1}}}+{\frac {\partial u_{2}}{\partial x_{1}}}{\frac {\partial u_{2}}{\partial x_{1}}}+{\frac {\partial u_{3}}{\partial x_{1}}}{\frac {\partial u_{3}}{\partial x_{1}}}\right)\\\varepsilon _{22}&={\frac {1}{2}}\left({\frac {\partial u_{2}}{\partial x_{2}}}+{\frac {\partial u_{2}}{\partial x_{2}}}\right)+{\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{2}}}{\frac {\partial u_{1}}{\partial x_{2}}}+{\frac {\partial u_{2}}{\partial x_{2}}}{\frac {\partial u_{2}}{\partial x_{2}}}+{\frac {\partial u_{3}}{\partial x_{2}}}{\frac {\partial u_{3}}{\partial x_{2}}}\right)\\\varepsilon _{33}&={\frac {1}{2}}\left({\frac {\partial u_{3}}{\partial x_{3}}}+{\frac {\partial u_{3}}{\partial x_{3}}}\right)+{\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{3}}}{\frac {\partial u_{1}}{\partial x_{3}}}+{\frac {\partial u_{2}}{\partial x_{3}}}{\frac {\partial u_{2}}{\partial x_{3}}}+{\frac {\partial u_{3}}{\partial x_{3}}}{\frac {\partial u_{3}}{\partial x_{3}}}\right)\\\varepsilon _{23}&={\frac {1}{2}}\left({\frac {\partial u_{2}}{\partial x_{3}}}+{\frac {\partial u_{3}}{\partial x_{2}}}\right)+{\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{2}}}{\frac {\partial u_{1}}{\partial x_{3}}}+{\frac {\partial u_{2}}{\partial x_{2}}}{\frac {\partial u_{2}}{\partial x_{3}}}+{\frac {\partial u_{3}}{\partial x_{2}}}{\frac {\partial u_{3}}{\partial x_{3}}}\right)\\\varepsilon _{31}&={\frac {1}{2}}\left({\frac {\partial u_{3}}{\partial x_{1}}}+{\frac {\partial u_{1}}{\partial x_{3}}}\right)+{\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{3}}}{\frac {\partial u_{1}}{\partial x_{1}}}+{\frac {\partial u_{2}}{\partial x_{3}}}{\frac {\partial u_{2}}{\partial x_{1}}}+{\frac {\partial u_{3}}{\partial x_{3}}}{\frac {\partial u_{3}}{\partial x_{1}}}\right)\\\varepsilon _{12}&={\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{2}}}+{\frac {\partial u_{2}}{\partial x_{1}}}\right)+{\frac {1}{2}}\left({\frac {\partial u_{1}}{\partial x_{1}}}{\frac {\partial u_{1}}{\partial x_{2}}}+{\frac {\partial u_{2}}{\partial x_{1}}}{\frac {\partial u_{2}}{\partial x_{2}}}+{\frac {\partial u_{3}}{\partial x_{1}}}{\frac {\partial u_{3}}{\partial x_{2}}}\right)\end{aligned}} The displacements

$u_{1}=u_{0}(x_{1})-x_{3}{\cfrac {dw_{0}}{dx_{1}}}~;~~u_{2}=0~;~~u_{3}=w_{0}(x_{1})$ The derivatives

{\begin{aligned}{\frac {\partial u_{1}}{\partial x_{1}}}&={\cfrac {du_{0}}{dx_{1}}}-x_{3}{\cfrac {d^{2}w_{0}}{dx_{1}^{2}}}~;~~&{\frac {\partial u_{1}}{\partial x_{2}}}&=0~;~~&{\frac {\partial u_{1}}{\partial x_{3}}}&=-{\cfrac {dw_{0}}{dx_{1}}}\\{\frac {\partial u_{2}}{\partial x_{1}}}&=0~;~~&{\frac {\partial u_{2}}{\partial x_{2}}}&=0~;~~&{\frac {\partial u_{2}}{\partial x_{3}}}&=0\\{\frac {\partial u_{3}}{\partial x_{1}}}&={\cfrac {dw_{0}}{dx_{1}}}~;~~&{\frac {\partial u_{3}}{\partial x_{2}}}&=0~;~~&{\frac {\partial u_{3}}{\partial x_{3}}}&=0\end{aligned}} ### von Karman strains

The von Karman strains

{\begin{aligned}\varepsilon _{11}&={\cfrac {du_{0}}{dx_{1}}}-x_{3}{\cfrac {d^{2}w_{0}}{dx_{1}^{2}}}+{\frac {1}{2}}\left[\left({\cfrac {du_{0}}{dx_{1}}}-x_{3}{\cfrac {d^{2}w_{0}}{dx_{1}^{2}}}\right)^{2}+\left({\cfrac {dw_{0}}{dx_{1}}}\right)^{2}\right]\\\varepsilon _{22}&=0\\\varepsilon _{33}&={\frac {1}{2}}\left({\cfrac {dw_{0}}{dx_{1}}}\right)^{2}\\\varepsilon _{23}&=0\\\varepsilon _{31}&={\frac {1}{2}}\left({\cfrac {dw_{0}}{dx_{1}}}-{\cfrac {dw_{0}}{dx_{1}}}\right)-{\frac {1}{2}}\left[\left({\cfrac {du_{0}}{dx_{1}}}-x_{3}{\cfrac {d^{2}w_{0}}{dx_{1}^{2}}}\right)\left({\cfrac {dw_{0}}{dx_{1}}}\right)\right]\\\varepsilon _{12}&=0\end{aligned}} ## Equilibrium Equations

### Balance of forces

{\begin{aligned}{\cfrac {dN_{xx}}{dx}}+f(x)&=0\\{\cfrac {d^{2}M_{xx}}{dx^{2}}}+q(x)+{\cfrac {d}{dx}}\left(N_{xx}{\cfrac {dw_{0}}{dx}}\right)&=0\end{aligned}} ### Stress Resultants

{\begin{aligned}N_{xx}&=\int _{A}\sigma _{xx}~dA\\M_{xx}&=\int _{A}z\sigma _{xx}~dA\end{aligned}} ## Constitutive Relations

### Stress-Strain equation

$\sigma _{xx}=E\varepsilon _{xx}$ ### Stress Resultant - Displacement relations

{\begin{aligned}N_{xx}&=A_{xx}\varepsilon _{xx}^{0}+B_{xx}\varepsilon _{xx}^{1}=A_{xx}\left[{\cfrac {du_{0}}{dx}}+{\frac {1}{2}}\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]-B_{xx}{\cfrac {d^{2}w_{0}}{dx^{2}}}\\M_{xx}&=B_{xx}\varepsilon _{xx}^{0}+D_{xx}\varepsilon _{xx}^{1}=B_{xx}\left[{\cfrac {du_{0}}{dx}}+{\frac {1}{2}}\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]-D_{xx}{\cfrac {d^{2}w_{0}}{dx^{2}}}\end{aligned}} ### Extensional/Bending Stiffness

{\begin{aligned}A_{xx}&=\int _{A}E~dA\qquad \leftarrow \qquad {\text{extensional stiffness}}\\B_{xx}&=\int _{A}zE~dA\qquad \leftarrow \qquad {\text{extensional-bending stiffness}}\\D_{xx}&=\int _{A}z^{2}E~dA\qquad \leftarrow \qquad {\text{bending stiffness}}\end{aligned}} If $E$ is constant, and $x$ -axis passes through centroid

{\begin{aligned}A_{xx}&=E\int _{A}~dA=EA\\B_{xx}&=E\int _{A}z~dA=0\\D_{xx}&=E\int _{A}z^{2}~dA=EI\end{aligned}} ## Weak Forms

### Axial Equation

{\begin{aligned}\int _{x_{a}}^{x_{b}}{\cfrac {d(\delta u_{0})}{dx}}\left[{\cfrac {du_{0}}{dx}}+{\frac {1}{2}}\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]A_{xx}~dx&=\int _{x_{a}}^{x_{b}}(\delta u_{0})f~dx+\\&\delta u_{0}(x_{a})Q_{1}+\delta u_{0}(x_{b})Q_{4}\end{aligned}} where

{\begin{aligned}\delta u_{0}&:=v_{1}\\Q_{1}&:=-N_{xx}(x_{a})\\Q_{4}&:=N_{xx}(x_{b})\end{aligned}} ### Bending Equation

{\begin{aligned}\int _{x_{a}}^{x_{b}}\left\{{\cfrac {d(\delta w_{0})}{dx}}\right.&\left[{\cfrac {du_{0}}{dx}}+{\cfrac {1}{2}}~\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]{\cfrac {dw_{0}}{dx}}A_{xx}+\left.{\cfrac {d^{2}(\delta w_{0})}{dx^{2}}}\left({\cfrac {d^{2}w_{0}}{dx^{2}}}\right)D_{xx}\right\}~dx=\\&\int _{x_{a}}^{x_{b}}(\delta w_{0})q~dx+\delta w_{0}(x_{a})Q_{2}+\delta w_{0}(x_{b})Q_{5}+\delta \theta (x_{a})Q_{3}+\delta \theta (x_{b})Q_{6}~.\end{aligned}} where

{\begin{aligned}\delta w_{0}&:=v_{2}&\delta \theta &:={\cfrac {dv_{2}}{dx}}\\Q_{2}&:=-\left[{\cfrac {dM_{xx}}{dx}}+N_{xx}{\cfrac {dw_{0}}{dx}}\right]_{x_{a}}&Q_{5}&:=\left[{\cfrac {dM_{xx}}{dx}}+N_{xx}{\cfrac {dw_{0}}{dx}}\right]_{x_{b}}\\Q_{3}&:=-M_{xx}(x_{a})&Q_{6}&:=M_{xx}(x_{b})\end{aligned}} ## Finite Element Model Finite element model for Euler Bernoulli beam
{\begin{aligned}u_{0}(x)&=u_{1}\psi _{1}(x)+u_{2}\psi _{2}(x)\\w_{0}(x)&=w_{1}\phi _{1}(x)+\theta _{1}\phi _{2}(x)+w_{2}\phi _{3}(x)+\theta _{2}\phi _{4}(x)\end{aligned}} where $\theta =-(dw_{0}/dx)$ .

### Hermite Cubic Shape Functions Hermite shape functions for beam

### Finite Element Equations

${\begin{bmatrix}\mathbf {K} ^{11}&\vdots &\mathbf {K} ^{12}\\&\vdots &\\\mathbf {K} ^{21}&\vdots &\mathbf {K} ^{22}\end{bmatrix}}{\begin{bmatrix}\mathbf {u} \\\\\mathbf {d} \end{bmatrix}}={\begin{bmatrix}\mathbf {F} ^{1}\\\\\mathbf {F} ^{2}\end{bmatrix}}$ where

{\begin{aligned}\mathbf {u} &=[u_{1}\quad u_{2}]^{T}\\\mathbf {d} &=[w_{1}\quad \theta _{1}\quad w_{2}\quad \theta _{2}]^{T}\end{aligned}} {\begin{aligned}\mathbf {K} ^{11}&=2\times 2;\qquad &\mathbf {K} ^{12}&=2\times 4\\\mathbf {K} ^{21}&=4\times 2;\qquad &\mathbf {K} ^{22}&=4\times 4\end{aligned}} ### Symmetric Stiffness Matrix

{\begin{aligned}K_{ij}^{11}&=\int _{x_{a}}^{x_{b}}A_{xx}{\cfrac {d\psi _{i}}{dx}}{\cfrac {d\psi _{j}}{dx}}~dx\\K_{ij}^{12}&={\frac {1}{2}}\int _{x_{a}}^{x_{b}}\left(A_{xx}{\cfrac {dw_{0}}{dx}}\right){\cfrac {d\psi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}~dx\\K_{ij}^{21}&={\frac {1}{2}}\int _{x_{a}}^{x_{b}}\left(A_{xx}{\cfrac {dw_{0}}{dx}}\right){\cfrac {d\phi _{i}}{dx}}{\cfrac {d\psi _{j}}{dx}}~dx\\K_{ij}^{22}&=\int _{x_{a}}^{x_{b}}\left\{{\frac {1}{2}}A_{xx}\left[\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]{\cfrac {d\phi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}+D_{xx}{\cfrac {d^{2}\phi _{i}}{dx^{2}}}{\cfrac {d^{2}\phi _{j}}{dx^{2}}}\right\}~dx\end{aligned}} {\begin{aligned}F_{i}^{1}&=\int _{x_{a}}^{x_{b}}\psi _{i}f~dx+\psi _{i}(x_{a})Q_{1}+\psi _{i}(x_{b})Q_{4}\\F_{i}^{2}&=\int _{x_{a}}^{x_{b}}\phi _{i}q~dx+\phi _{i}(x_{a})Q_{2}+\phi _{i}(x_{b})Q_{5}+{\cfrac {d\phi _{i}}{dx}}(x_{a})Q_{3}+{\cfrac {d\phi _{i}}{dx}}(x_{b})Q_{6}\end{aligned}} ### Newton-Raphson Solution

$\mathbf {K} (\mathbf {U} )\mathbf {U} =\mathbf {F}$ where

{\begin{aligned}U_{1}&=u_{1},~U_{2}=u_{2},~U_{3}=d_{1},~U_{4}=d_{2},~U_{5}=d_{3},~U_{6}=d_{4}\\F_{1}&=F_{1}^{1},~F_{2}=F_{2}^{1},~F_{3}=F_{1}^{2},~F_{4}=F_{2}^{2},~F_{5}=F_{3}^{2},~F_{6}=F_{4}^{2}\end{aligned}} The residual is

$\mathbf {R} =\mathbf {K} \mathbf {U} -\mathbf {F} ~.$ For Newton iterations, we use the algorithm

$\mathbf {U} ^{r+1}=\mathbf {U} ^{r}-(\mathbf {T} ^{r})^{-1}\mathbf {R} ^{r}$ where the tangent stiffness matrix is given by

$\mathbf {T} ^{r}={\frac {\partial \mathbf {R} ^{r}}{\partial \mathbf {U} }};\quad {\text{or}}\quad T_{ij}={\frac {\partial R_{i}}{\partial U_{j}}},\qquad i=1\dots 6,j=1\dots 6~.$ ### Tangent Stiffness Matrix

{\begin{aligned}i=1\dots 2;~j=1\dots 2&:\\&{T_{ij}^{11}=K_{ij}^{11}}\\\\i=1\dots 2;~j=1\dots 4&:\\&{T_{ij}^{12}=2K_{ij}^{12}}\\\\i=1\dots 4;~j=1\dots 2&:\\&{T_{ij}^{21}=2K_{ij}^{21}}\\\\i=1\dots 4;~j=1\dots 4&:\\&{T_{ij}^{22}=K_{ij}^{22}+{\frac {1}{2}}\int _{x_{a}}^{x_{b}}A_{xx}\left[{\cfrac {du_{0}}{dx}}+2\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]{\cfrac {d\phi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}~dx}\end{aligned}} Recall

$N_{xx}=A_{xx}\left[{\cfrac {du_{0}}{dx}}+{\frac {1}{2}}\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]-{B_{xx}}~~{\cfrac {d^{2}w_{0}}{dx^{2}}}$ • Divide load into small increments.
$F=\sum _{i=1}^{N}\Delta F_{i}$ • Compute $\mathbf {u}$ and $\mathbf {d}$ for first load step,
$\mathbf {K} (\mathbf {U} _{0})\mathbf {U} _{1}=\Delta F_{1}$  Stiffness of Euler-Bernoulli beam.
• Compute $\mathbf {u}$ and $\mathbf {d}$ for second load step,
$\mathbf {K} (\mathbf {U} _{1})\mathbf {U} _{2}=\Delta F_{1}+\Delta F_{2}$ • Continue until F is reached.

## Membrane Locking

Recall

${\begin{bmatrix}\mathbf {K} ^{11}&\vdots &\mathbf {K} ^{12}\\&\vdots &\\\mathbf {K} ^{21}&\vdots &\mathbf {K} ^{22}\end{bmatrix}}{\begin{bmatrix}\mathbf {u} \\\\\mathbf {d} \end{bmatrix}}={\begin{bmatrix}\mathbf {F} ^{1}\\\\\mathbf {F} ^{2}\end{bmatrix}}$ where

{\begin{aligned}K_{ij}^{11}&=\int _{x_{a}}^{x_{b}}A_{xx}{\cfrac {d\psi _{i}}{dx}}{\cfrac {d\psi _{j}}{dx}}~dx\\K_{ij}^{12}&={\frac {1}{2}}\int _{x_{a}}^{x_{b}}\left(A_{xx}{\cfrac {dw_{0}}{dx}}\right){\cfrac {d\psi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}~dx\\K_{ij}^{21}&={\frac {1}{2}}\int _{x_{a}}^{x_{b}}\left(A_{xx}{\cfrac {dw_{0}}{dx}}\right){\cfrac {d\phi _{i}}{dx}}{\cfrac {d\psi _{j}}{dx}}~dx\\K_{ij}^{22}&=\int _{x_{a}}^{x_{b}}\left\{{\frac {1}{2}}A_{xx}\left[\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]{\cfrac {d\phi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}+D_{xx}{\cfrac {d^{2}\phi _{i}}{dx^{2}}}{\cfrac {d^{2}\phi _{j}}{dx^{2}}}\right\}~dx\end{aligned}}  Mebrane locking in Euler-Bernoulli beam

### For Hinged-Hinged

Membrane strain:

$\varepsilon _{xx}^{0}=0$ or

${\cfrac {du_{0}}{dx}}+{\frac {1}{2}}\left({\cfrac {dw_{0}}{dx}}\right)^{2}=0$ Hence, shape functions should be such that

${\cfrac {du_{0}}{dx}}\approx \left({\cfrac {dw_{0}}{dx}}\right)^{2}$ $u_{0}$ linear, $w_{0}$ cubic $\implies$ Element Locks! Too stiff.

## Selective Reduced Integration

• Assume $u_{0}$ is linear ;~~ $w_{0}$ is cubic.
• Then ${\cfrac {du_{0}}{dx}}\equiv {\cfrac {d\psi _{i}}{dx}}$ is constant, and ${\cfrac {dw_{0}}{dx}}\equiv {\cfrac {d\phi _{i}}{dx}}$ is quadratic.
• Try to keep $\varepsilon _{xx}^{0}=$ constant.
{\begin{aligned}K_{ij}^{11}&=\int _{x_{a}}^{x_{b}}A_{xx}{\cfrac {d\psi _{i}}{dx}}{\cfrac {d\psi _{j}}{dx}}~dx\\K_{ij}^{12}&={\frac {1}{2}}\int _{x_{a}}^{x_{b}}\left(A_{xx}{\cfrac {dw_{0}}{dx}}\right){\cfrac {d\psi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}~dx\\K_{ij}^{22}&=\int _{x_{a}}^{x_{b}}\left\{{\frac {1}{2}}A_{xx}\left[{\cfrac {du_{0}}{dx}}+\left({\cfrac {dw_{0}}{dx}}\right)^{2}\right]{\cfrac {d\phi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}+D_{xx}{\cfrac {d^{2}\phi _{i}}{dx^{2}}}{\cfrac {d^{2}\phi _{j}}{dx^{2}}}\right\}~dx\end{aligned}} • $\mathbf {K} ^{11}$ integrand is constant, $\mathbf {K} ^{12}$ integrand is fourth-order , $\mathbf {K} ^{22}$ integrand is eighth-order

### Full integration

$n_{\text{gauss pt}}={\text{int}}[(p+1)/2]+1$ Assume $A_{xx}$ = constant.

{\begin{aligned}K_{ij}^{11}&=A_{xx}\int _{x_{a}}^{x_{b}}{\cfrac {d\psi _{i}}{dx}}{\cfrac {d\psi _{j}}{dx}}~dx\\&=A_{xx}\int _{-1}^{1}\left(J^{-1}{\cfrac {d\psi _{i}(\xi )}{d\xi }}\right)\left(J^{-1}{\cfrac {d\psi _{j}(\xi )}{d\xi }}\right)J~d\xi =A_{xx}\int _{-1}^{1}F(\xi )~d\xi \\&\approx A_{xx}W_{1}F(\xi _{1})\leftarrow {\text{one-point integration}}\end{aligned}} {\begin{aligned}K_{ij}^{12}&={\frac {1}{2}}\int _{x_{a}}^{x_{b}}\left(A_{xx}{\cfrac {dw_{0}}{dx}}\right){\cfrac {d\psi _{i}}{dx}}{\cfrac {d\phi _{j}}{dx}}~dx\\&={\cfrac {A_{xx}}{2}}\int _{-1}^{1}\left(\sum _{i=1}^{4}w_{i}J^{-1}{\cfrac {d\phi _{i}(\xi )}{d\xi }}\right)\left(J^{-1}{\cfrac {d\psi _{i}(\xi )}{d\xi }}\right)\left(J^{-1}{\cfrac {d\phi _{j}(\xi )}{d\xi }}\right)~dx\\&\approx A_{xx}\left[W_{1}F(\xi _{1})+W_{2}F(\xi _{2})+W_{3}F(\xi _{3})\right]\leftarrow {\text{full integration}}\\&\approx A_{xx}\left[W_{1}F(\xi _{1})+W_{2}F(\xi _{2})\right]\leftarrow {\text{reduced integration}}\end{aligned}} 