Nonlinear finite elements/Homework 5/Solutions/Problem 1

Problem 1: Nonlinear Beam Bending

Given:

The differential equations governing the bending of straight beams are

$\begin{array}{cc}\frac{d{N}_{xx}}{dx}+f(x)& =0\\ \frac{dV}{dx}+q(x)& =0\\ \frac{d{M}_{xx}}{dx}-V+N_{xx}\frac{d{w}_0}{dx}& =0\end{array}$

Solution

Part 1

Show that the weak forms of these equations can be written as

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_1}{dx}{N}_{xx}dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{1}fdx+{v_1{N}_{xx}|}_{x_a}^{x_b}\\ \underset{x_a}{\overset{x_b}{\int }}[\frac{d{v}_2}{dx}\left(\frac{d{w}_0}{dx}{N}_{xx}\right)-\frac{d^2{v}_2}{d{x}^2}{M}_{xx}]dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{2}qdx+{v_2(\frac{d{w}_0}{dx}{N}_{xx}+\frac{d{M}_{xx}}{dx})|}_{x_a}^{x_b}-{\frac{d{v}_2}{dx}{M}_{xx}|}_{x_a}^{x_b}\end{array}$

First we get rid of the shear force term by combining the second and third equations to get

$\begin{array}{cc}\frac{d{N}_{xx}}{dx}+f(x)& =0\text{(1)}\\ \frac{d^2{M}_{xx}}{d{x}^2}+q(x)+\frac{d}{dx}\left(N_{xx}\frac{d{w}_0}{dx}\right)& =0\text{(2)}\end{array}$

Let $v_1$ be the weighting function for equation (1) and let $v_2$ be the weighting function for equation (2).

Then the weak forms of the two equations are

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}{v}_1(\frac{d{N}_{xx}}{dx}+f(x))dx& =0\text{(3)}\\ \underset{x_a}{\overset{x_b}{\int }}{v}_2[\frac{d^2{M}_{xx}}{d{x}^2}+q(x)+\frac{d}{dx}\left(N_{xx}\frac{d{w}_0}{dx}\right)]dx& =0\text{(4)}\end{array}$

To get the symmetric weak forms, we integrate by parts (even though the symmetry is not obvious at this stage) to get

$\begin{array}{cc}{(v_1{N}_{xx})|}_{x_a}^{x_b}& -\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_1}{dx}{N}_{xx}dx+\underset{x_a}{\overset{x_b}{\int }}{v}_{1}f(x)dx=0\text{(5)}\\ {\left(v_2\frac{d{M}_{xx}}{dx}\right)|}_{x_a}^{x_b}& -\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_2}{dx}\frac{d{M}_{xx}}{dx}dx+\underset{x_a}{\overset{x_b}{\int }}{v}_{2}q(x)dx+\\ & {\left(v_2{N}_{xx}\frac{d{w}_0}{dx}\right)|}_{x_a}^{x_b}-\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_2}{dx}\left(N_{xx}\frac{d{w}_0}{dx}\right)dx=0\text{(6)}\end{array}$

Integrate equation (6) again by parts, and get

$\begin{array}{cc}{\left(v_2\frac{d{M}_{xx}}{dx}\right)|}_{x_a}^{x_b}& -{\left(\frac{d{v}_2}{dx}{M}_{xx}\right)|}_{x_a}^{x_b}+\underset{x_a}{\overset{x_b}{\int }}\frac{d^2{v}_2}{d{x}^2}{M}_{xx}dx+\underset{x_a}{\overset{x_b}{\int }}{v}_{2}q(x)dx+\\ & {\left(v_2{N}_{xx}\frac{d{w}_0}{dx}\right)|}_{x_a}^{x_b}-\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_2}{dx}\left(N_{xx}\frac{d{w}_0}{dx}\right)dx=0\text{(7)}\end{array}$

Collect terms and rearrange equations (5) and (7) to get

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_1}{dx}{N}_{xx}dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{1}f(x)dx+{(v_1{N}_{xx})|}_{x_a}^{x_b}\text{(8)}\\ \underset{x_a}{\overset{x_b}{\int }}[\frac{d{v}_2}{dx}\left(N_{xx}\frac{d{w}_0}{dx}\right)-\frac{d^2{v}_2}{d{x}^2}{M}_{xx}]dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{2}q(x)dx+\\ & {[\left(v_2\frac{d{M}_{xx}}{dx}\right)-\left(\frac{d{v}_2}{dx}{M}_{xx}\right)+\left(v_2{N}_{xx}\frac{d{w}_0}{dx}\right)]}_{x_a}^{x_b}\text{(9)}\end{array}$

Rewriting the equations, we get

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_1}{dx}{N}_{xx}dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{1}fdx+{v_1{N}_{xx}|}_{x_a}^{x_b}\\ \underset{x_a}{\overset{x_b}{\int }}[\frac{d{v}_2}{dx}\left(N_{xx}\frac{d{w}_0}{dx}\right)-\frac{d^2{v}_2}{d{x}^2}{M}_{xx}]dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{2}qdx+{\left[v_2(\frac{d{M}_{xx}}{dx}+N_{xx}\frac{d{w}_0}{dx})\right]}_{x_a}^{x_b}-{\left[\frac{d{v}_2}{dx}{M}_{xx}\right]}_{x_a}^{x_b}\end{array}\text{(10)}$

Hence shown.

Part 2

The von Karman strains are related to the displacements by

$\begin{array}{cc}{\epsilon }_{xx}& ={\epsilon }_{xx}^0+z{\epsilon }_{xx}^1\\ {\epsilon }_{xx}^0& =\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2\\ {\epsilon }_{xx}^1& =-\frac{d^2{w}_0}{d{x}^2}\end{array}$

The stress and moment resultants are defined as

$\begin{array}{cc}{N}_{xx}& =\underset{A}{\int }{\sigma }_{xx}dA\\ M_{xx}& =\underset{A}{\int }z{\sigma }_{xx}dA\end{array}$

For a linear elastic material, the stiffnesses of the beam in extension and bending are defined as

$\begin{array}{cc}{A}_{xx}& =\underset{A}{\int }EdA\leftarrow \text{extensional stiffness}\\ B_{xx}& =\underset{A}{\int }zEdA\leftarrow \text{extensional-bending stiffness}\\ D_{xx}& =\underset{A}{\int }{z}^{2}EdA\leftarrow \text{bending stiffness}\end{array}$

where $E$ is the Young's modulus of the material.

Derive expressions for $N_{xx}$ and $M_{xx}$ in terms of the displacements $u_0$ and $w_0$ and the extensional and bending stiffnesses of the beam assuming a linear elastic material.

The stress-strain relations for an isotropic linear elastic material are

$\left[\begin{array}{c}{\sigma }_{xx}\\ {\sigma }_{yy}\\ {\sigma }_{zz}\\ {\sigma }_{yz}\\ {\sigma }_{zx}\\ {\sigma }_{xy}\end{array}\right]=\frac{E}{(1+\nu )(1-2\nu )}\left[\begin{array}{cccccc}1-\nu & \nu & \nu & 0& 0& 0\\ \nu & 1-\nu & \nu & 0& 0& 0\\ \nu & \nu & 1-\nu & 0& 0& 0\\ 0& 0& 0& (1-2\nu )& 0& 0\\ 0& 0& 0& 0& (1-2\nu )& 0\\ 0& 0& 0& 0& 0& (1-2\nu )\end{array}\right]\left[\begin{array}{c}{\epsilon }_{xx}\\ {\epsilon }_{yy}\\ {\epsilon }_{zz}\\ {\epsilon }_{yz}\\ {\epsilon }_{zx}\\ {\epsilon }_{xy}\end{array}\right]$

Since all strains other than ${\epsilon }_{xx}$ are zero, the above equations reduce to

${\sigma }_{xx}=\frac{E(1-\nu )}{(1+\nu )(1-2\nu )}{\epsilon }_{xx};{\sigma }_{yy}=\frac{E\nu }{(1+\nu )(1-2\nu )}{\epsilon }_{xx};{\sigma }_{zz}=\frac{E\nu }{(1+\nu )(1-2\nu )}{\epsilon }_{xx}.$

If we ignore the stresses ${\sigma }_{yy}$ and ${\sigma }_{zz}$, the only allowable value of $\nu$ is zero. Then the stress-strain relations become

${\sigma }_{xx}=E{\epsilon }_{xx}.$

Plugging this relation into the stress and moment of stress resultant equations, we get

$\begin{array}{cc}{N}_{xx}& =\underset{A}{\int }E{\epsilon }_{xx}dA\\ M_{xx}& =\underset{A}{\int }zE{\epsilon }_{xx}dA\end{array}$

Plugging in the relations for the strain we get

$\begin{array}{cc}{N}_{xx}& =\underset{A}{\int }E({\epsilon }_{xx}^0+z{\epsilon }_{xx}^1)dA=\underset{A}{\int }E{\epsilon }_{xx}^{0}dA+\underset{A}{\int }zE{\epsilon }_{xx}^{1}dA\\ M_{xx}& =\underset{A}{\int }zE({\epsilon }_{xx}^0+z{\epsilon }_{xx}^1)dA=\underset{A}{\int }zE{\epsilon }_{xx}^{0}dA+\underset{A}{\int }{z}^{2}E{\epsilon }_{xx}^{1}dA.\end{array}$

Since both ${\epsilon }_{xx}^0$ and ${\epsilon }_{xx}^1$ are independent of $y$ and $z$, we can take these quantities outside the integrals and get

$\begin{array}{cc}{N}_{xx}& ={\epsilon }_{xx}^0\underset{A}{\int }EdA+{\epsilon }_{xx}^1\underset{A}{\int }zEdA\\ M_{xx}& ={\epsilon }_{xx}^0\underset{A}{\int }zEdA+{\epsilon }_{xx}^1\underset{A}{\int }{z}^{2}EdA.\end{array}$

Using the definitions of the extensional, extensional-bending, and bending stiffness, we can then write

$\begin{array}{cc}{N}_{xx}& ={\epsilon }_{xx}^0{A}_{xx}+{\epsilon }_{xx}^1{B}_{xx}\\ M_{xx}& ={\epsilon }_{xx}^0{B}_{xx}+{\epsilon }_{xx}^1{D}_{xx}.\end{array}$

To write these relations in terms of $u_0$ and $w_0$ we substitute the expressions for the von Karman strains to get

$\begin{array}{cc}{N}_{xx}& =[\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]A_{xx}-\frac{d^2{w}_0}{d{x}^2}{B}_{xx}\\ M_{xx}& =[\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]B_{xx}-\frac{d^2{w}_0}{d{x}^2}{D}_{xx}.\end{array}$

These are the expressions of the resultants in terms of the displacements.

Part 3

Express the weak forms in terms of the displacements and the extensional and bending stiffnesses.

The weak form equations are

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_1}{dx}{N}_{xx}dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{1}fdx+{v_1{N}_{xx}|}_{x_a}^{x_b}\text{(11)}\\ \underset{x_a}{\overset{x_b}{\int }}[\frac{d{v}_2}{dx}\left(N_{xx}\frac{d{w}_0}{dx}\right)-\frac{d^2{v}_2}{d{x}^2}{M}_{xx}]dx& =\underset{x_a}{\overset{x_b}{\int }}{v}_{2}qdx+{\left[v_2(\frac{d{M}_{xx}}{dx}+N_{xx}\frac{d{w}_0}{dx})\right]}_{x_a}^{x_b}-{\left[\frac{d{v}_2}{dx}{M}_{xx}\right]}_{x_a}^{x_b}\text{(12)}\end{array}$

At this stage we make two more assumptions:

1. The elastic modulus is constant throughout the cross-section.
2. The $x$-axis passes through the centroid of the cross-section.

From the first assumption, we have

$B_{xx}=E\underset{A}{\int }zdA.$

From the second assumption, we get

$\underset{A}{\int }zdA=0⟹B_{xx}=0.$

Then the relations for $N_{xx}$ and $M_{xx}$ reduce to

$\begin{array}{cc}{N}_{xx}& =[\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]A_{xx}\text{(13)}\\ M_{xx}& =-\frac{d^2{w}_0}{d{x}^2}{D}_{xx}\text{(14)}.\end{array}$

Let us first consider equation (11). Plugging in the expression for $N_{xx}$ we get

$\underset{x_a}{\overset{x_b}{\int }}\frac{d{v}_1}{dx}[\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]A_{xx}dx=\underset{x_a}{\overset{x_b}{\int }}{v}_{1}fdx+{v_1{N}_{xx}|}_{x_a}^{x_b}$

We can also write the above in terms of virtual displacements by defining $\delta u_0:=v_1$, $Q_1:=-N_{xx}(x_a)$, and $Q_4:=N_{xx}(x_b)$. Then we get

$\underset{x_a}{\overset{x_b}{\int }}\frac{d(\delta u_0)}{dx}[\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]A_{xx}dx=\underset{x_a}{\overset{x_b}{\int }}(\delta u_0)fdx+\delta u_0(x_a)Q_1+\delta u_0(x_b)Q_4.$

Next, we do the same for equation (12). Plugging in the expressions for $N_{xx}$ and $M_{xx}$, we get

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\{\frac{d{v}_2}{dx}\left([\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]A_{xx}\frac{d{w}_0}{dx}\right)& +\frac{d^2{v}_2}{d{x}^2}\left(\frac{d^2{w}_0}{d{x}^2}\right)D_{xx}\}dx=\underset{x_a}{\overset{x_b}{\int }}{v}_{2}qdx+\\ & {\left[v_2(\frac{d{M}_{xx}}{dx}+N_{xx}\frac{d{w}_0}{dx})\right]}_{x_a}^{x_b}-{\left[\frac{d{v}_2}{dx}{M}_{xx}\right]}_{x_a}^{x_b}\end{array}$

We can write the above in terms of the virtual displacements and the generalized forces by defining

$\begin{array}{cc}\delta w_0& :=v_2\\ \delta \theta & :=\frac{d{v}_2}{dx}\\ Q_2& :=-{[\frac{d{M}_{xx}}{dx}+N_{xx}\frac{d{w}_0}{dx}]}_{x_a}\\ Q_5& :={[\frac{d{M}_{xx}}{dx}+N_{xx}\frac{d{w}_0}{dx}]}_{x_b}\\ Q_3& :=-M_{xx}(x_a)\\ Q_6& :=M_{xx}(x_b)\end{array}$

to get

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\{\frac{d(\delta w_0)}{dx}& [\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]\frac{d{w}_0}{dx}{A}_{xx}+\frac{d^2(\delta w_0)}{d{x}^2}\left(\frac{d^2{w}_0}{d{x}^2}\right)D_{xx}\}dx=\\ & \underset{x_a}{\overset{x_b}{\int }}(\delta w_0)qdx+\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{array}$

Part 4

Assume that the approximate solutions for the axial displacement and the transverse deflection over a two noded element are given by

$\begin{array}{cc}{u}_0(x)& =u_1{\psi }_1(x)+u_2{\psi }_2(x)\text{(15)}\\ w_0(x)& =w_1{\varphi }_1(x)+{\theta }_1{\varphi }_2(x)+w_2{\varphi }_3(x)+{\theta }_2{\varphi }_4(x)\text{(16)}\end{array}$

where $\theta =-(d{w}_0/dx)$.

Compute the element stiffness matrix for the element.

The weak forms of the governing equations are

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\frac{d(\delta u_0)}{dx}& [\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]A_{xx}dx=\underset{x_a}{\overset{x_b}{\int }}(\delta u_0)fdx+\delta u_0(x_a)Q_1+\delta u_0(x_b)Q_4\text{(17)}\\ \underset{x_a}{\overset{x_b}{\int }}\{\frac{d(\delta w_0)}{dx}& [\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2]\frac{d{w}_0}{dx}{A}_{xx}+\frac{d^2(\delta w_0)}{d{x}^2}\left(\frac{d^2{w}_0}{d{x}^2}\right)D_{xx}\}dx=\\ & \underset{x_a}{\overset{x_b}{\int }}(\delta w_0)qdx+\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.\text{(18)}\end{array}$

Let us first write the approximate solutions as

$\begin{array}{cc}{u}_0(x)& =\sum _{j=1}^2{u}_j{\psi }_j\text{(19)}\\ w_0(x)& =\sum _{j=1}^4{d}_j{\varphi }_j\text{(20)}\end{array}$

where $d_j$ are generalized displacements and

$\{d_1,d_2,d_3,d_4\}\equiv \{w_1,{\theta }_1,w_2,{\theta }_2\}.$

To formulate the finite element system of equations, we substitute the expressions from $u_0$ and $w_0$ from equations (19) and (20) into the weak form, and substitute the shape functions ${\psi }_i$ for $\delta u_0$, ${\varphi }_i$ for $\delta w_0$.

For the first equation (17) we get

$\underset{x_a}{\overset{x_b}{\int }}\frac{d{\psi }_i}{dx}[\left(\sum _{j=1}^2{u}_j\frac{d{\psi }_j}{dx}\right)+\frac{1}{2}\frac{d{w}_0}{dx}\left(\sum _{j=1}^4{d}_j\frac{d{\varphi }_j}{dx}\right)]A_{xx}dx=\underset{x_a}{\overset{x_b}{\int }}{\psi }_{i}fdx+{\psi }_i(x_a)Q_1+{\psi }_i(x_b)Q_4.$

After reorganizing, we have

$\begin{array}{cc}\sum _{j=1}^2\left[\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\psi }_i}{dx}\frac{d{\psi }_j}{dx}dx\right]u_j& +\sum _{j=1}^4\left[\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\psi }_i}{dx}\frac{d{\varphi }_j}{dx}dx\right]d_j=\\ & \underset{x_a}{\overset{x_b}{\int }}{\psi }_{i}fdx+{\psi }_i(x_a)Q_1+{\psi }_i(x_b)Q_4.\end{array}$

We can write the above as

$\sum _{j=1}^2{K}_{ij}^{11}{u}_j+\sum _{j=1}^4{K}_{ij}^{12}{d}_j=F_i^1$

where $i=1,2$ and

$\begin{array}{cc}{K}_{ij}^{11}& =\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\psi }_i}{dx}\frac{d{\psi }_j}{dx}dx\\ K_{ij}^{12}& =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\psi }_i}{dx}\frac{d{\varphi }_j}{dx}dx\\ F_i^1& =\underset{x_a}{\overset{x_b}{\int }}{\psi }_{i}fdx+{\psi }_i(x_a)Q_1+{\psi }_i(x_b)Q_4.\end{array}$

For the second equation (18) we get

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\{\frac{d{\varphi }_i}{dx}& [\left(\sum _{j=1}^2{u}_j\frac{d{\psi }_j}{dx}\right)+\frac{1}{2}\frac{d{w}_0}{dx}\left(\sum _{j=1}^4{d}_j\frac{d{\varphi }_j}{dx}\right)]\frac{d{w}_0}{dx}{A}_{xx}+\frac{d^2{\varphi }_i}{d{x}^2}\left(\sum _{j=1}^4{d}_j\frac{d^2{\varphi }_j}{d{x}^2}\right)D_{xx}\}dx=\\ & \underset{x_a}{\overset{x_b}{\int }}{\varphi }_{i}qdx+{\varphi }_i(x_a)Q_2+{\varphi }_i(x_b)Q_5+\frac{d{\varphi }_i}{dx}(x_a)Q_3+\frac{d{\varphi }_i}{dx}(x_b)Q_6\end{array}$

After rearranging we get

$\begin{array}{cc}\sum _{j=1}^2\left[\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\varphi }_i}{dx}\frac{d{\psi }_j}{dx}dx\right]u_j& +\sum _{j=1}^4\left\{\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left[A_{xx}{\left(\frac{d{w}_0}{dx}\right)}^2\right]\frac{d{\varphi }_i}{dx}\frac{d{\varphi }_j}{dx}dx\right\}{d}_j\\ & +\sum _{j=1}^4\left(\underset{x_a}{\overset{x_b}{\int }}{D}_{xx}\frac{d^2{\varphi }_i}{d{x}^2}\frac{d^2{\varphi }_j}{d{x}^2}dx\right)d_j=\\ & \underset{x_a}{\overset{x_b}{\int }}{\varphi }_{i}qdx+{\varphi }_i(x_a)Q_2+{\varphi }_i(x_b)Q_5+\frac{d{\varphi }_i}{dx}(x_a)Q_3+\frac{d{\varphi }_i}{dx}(x_b)Q_6\end{array}$

We can write the above as

$\sum _{j=1}^2{K}_{ij}^{21}{u}_j+\sum _{j=1}^4{K}_{ij}^{22}{d}_j=F_i^2$

where $i=1,2,3,4$ and

$\begin{array}{cc}{K}_{ij}^{21}& =\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\varphi }_i}{dx}\frac{d{\psi }_j}{dx}dx\\ K_{ij}^{22}& =\underset{x_a}{\overset{x_b}{\int }}\{\frac{1}{2}\left[A_{xx}{\left(\frac{d{w}_0}{dx}\right)}^2\right]\frac{d{\varphi }_i}{dx}\frac{d{\varphi }_j}{dx}+D_{xx}\frac{d^2{\varphi }_i}{d{x}^2}\frac{d^2{\varphi }_j}{d{x}^2}\}dx\\ F_i^2& =\underset{x_a}{\overset{x_b}{\int }}{\varphi }_{i}qdx+{\varphi }_i(x_a)Q_2+{\varphi }_i(x_b)Q_5+\frac{d{\varphi }_i}{dx}(x_a)Q_3+\frac{d{\varphi }_i}{dx}(x_b)Q_6.\end{array}$

In matrix form, we can write

$𝐊=\left[\begin{array}{ccccccc}{K}_{11}^{11}& K_{12}^{11}& ⋮& K_{11}^{12}& K_{12}^{12}& K_{13}^{12}& K_{14}^{12}\\ K_{21}^{11}& K_{22}^{11}& ⋮& K_{21}^{12}& K_{22}^{12}& K_{23}^{12}& K_{24}^{12}\\ & & \dots & & & & \\ K_{11}^{21}& K_{12}^{21}& ⋮& K_{11}^{22}& K_{12}^{22}& K_{13}^{12}& K_{14}^{12}\\ K_{21}^{21}& K_{22}^{21}& ⋮& K_{21}^{22}& K_{22}^{22}& K_{23}^{12}& K_{24}^{12}\\ K_{31}^{21}& K_{32}^{21}& ⋮& K_{31}^{22}& K_{32}^{22}& K_{33}^{22}& K_{34}^{22}\\ K_{41}^{21}& K_{42}^{21}& ⋮& K_{41}^{22}& K_{42}^{22}& K_{43}^{22}& K_{44}^{22}\end{array}\right]$

or

$𝐊=\left[\begin{array}{ccc}{𝐊}^{11}& ⋮& {𝐊}^{12}\\ & ⋮& \\ {𝐊}^{21}& ⋮& {𝐊}^{22}\end{array}\right].$

The finite element system of equations can then be written as

$\text{(21)}\left[\begin{array}{ccccccc}{K}_{11}^{11}& K_{12}^{11}& ⋮& K_{11}^{12}& K_{12}^{12}& K_{13}^{12}& K_{14}^{12}\\ K_{21}^{11}& K_{22}^{11}& ⋮& K_{21}^{12}& K_{22}^{12}& K_{23}^{12}& K_{24}^{12}\\ & & \dots & & & & \\ K_{11}^{21}& K_{12}^{21}& ⋮& K_{11}^{22}& K_{12}^{22}& K_{13}^{12}& K_{14}^{12}\\ K_{21}^{21}& K_{22}^{21}& ⋮& K_{21}^{22}& K_{22}^{22}& K_{23}^{12}& K_{24}^{12}\\ K_{31}^{21}& K_{32}^{21}& ⋮& K_{31}^{22}& K_{32}^{22}& K_{33}^{22}& K_{34}^{22}\\ K_{41}^{21}& K_{42}^{21}& ⋮& K_{41}^{22}& K_{42}^{22}& K_{43}^{22}& K_{44}^{22}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ \\ d_1\\ d_2\\ d_3\\ d_4\end{array}\right]=\left[\begin{array}{c}{F}_1^1\\ F_2^1\\ \\ F_1^2\\ F_2^2\\ F_3^2\\ F_4^2\end{array}\right]$

or

$\left[\begin{array}{ccc}{𝐊}^{11}& ⋮& {𝐊}^{12}\\ & ⋮& \\ {𝐊}^{21}& ⋮& {𝐊}^{22}\end{array}\right]\left[\begin{array}{c}𝐮\\ \\ 𝐝\end{array}\right]=\left[\begin{array}{c}{𝐅}^1\\ \\ {𝐅}^2\end{array}\right].$

Part 5

Show the alternate procedure by which the element stiffness matrix can be made symmetric.

The stiffness matrix is unsymmetric because ${𝐊}^{12}$ contains a factor of $1/2$ while ${𝐊}^{21}$ does not. The expressions of these terms are

$\begin{array}{cc}{K}_{ij}^{12}& =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\psi }_i}{dx}\frac{d{\varphi }_j}{dx}dx,i=1,2,\text{and}j=1,2,3,4\\ K_{ij}^{21}& =\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\varphi }_i}{dx}\frac{d{\psi }_j}{dx}dx,i=1,2,3,4,\text{and}j=1,2.\end{array}$

To get a symmetric stiffness matrix, we write equation (18) as

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\{\frac{d(\delta w_0)}{dx}& [\frac{1}{2}\frac{d{u}_0}{dx}+\frac{1}{2}{\left(\frac{d{w}_0}{dx}\right)}^2+\frac{1}{2}\frac{d{u}_0}{dx}]\frac{d{w}_0}{dx}{A}_{xx}+\frac{d^2(\delta w_0)}{d{x}^2}\left(\frac{d^2{w}_0}{d{x}^2}\right)D_{xx}\}dx=\\ & \underset{x_a}{\overset{x_b}{\int }}(\delta w_0)qdx+\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{array}$

The quantity is green is assumed to be known from a previous iteration and adds to the ${𝐊}^{22}$ terms.

Repeating the procedure used in the previous question

$\begin{array}{cc}\underset{x_a}{\overset{x_b}{\int }}\{\frac{d{\varphi }_i}{dx}& [\frac{1}{2}\left(\sum _{j=1}^2{u}_j\frac{d{\psi }_j}{dx}\right)+\frac{1}{2}\frac{d{w}_0}{dx}\left(\sum _{j=1}^4{d}_j\frac{d{\varphi }_j}{dx}\right)]\frac{d{w}_0}{dx}{A}_{xx}+\frac{1}{2}\frac{d{\varphi }_i}{dx}\frac{d{u}_0}{dx}\left(\sum _{j=1}^4{d}_j\frac{d{\varphi }_j}{dx}\right)A_{xx}+\\ & \frac{d^2{\varphi }_i}{d{x}^2}\left(\sum _{j=1}^4{d}_j\frac{d^2{\varphi }_j}{d{x}^2}\right)D_{xx}\}dx=\\ & \underset{x_a}{\overset{x_b}{\int }}{\varphi }_{i}qdx+{\varphi }_i(x_a)Q_2+{\varphi }_i(x_b)Q_5+\frac{d{\varphi }_i}{dx}(x_a)Q_3+\frac{d{\varphi }_i}{dx}(x_b)Q_6\end{array}$

After rearranging we get

$\begin{array}{cc}\sum _{j=1}^2& \left[\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\varphi }_i}{dx}\frac{d{\psi }_j}{dx}dx\right]u_j+\sum _{j=1}^4\left\{\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left[A_{xx}{\left(\frac{d{w}_0}{dx}\right)}^2\right]\frac{d{\varphi }_i}{dx}\frac{d{\varphi }_j}{dx}dx\right\}{d}_j+\\ & \sum _{j=1}^4\left[\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{u}_0}{dx}\frac{d{\varphi }_i}{dx}\frac{d{\varphi }_j}{dx}dx\right]d_j+\sum _{j=1}^4\left(\underset{x_a}{\overset{x_b}{\int }}{D}_{xx}\frac{d^2{\varphi }_i}{d{x}^2}\frac{d^2{\varphi }_j}{d{x}^2}dx\right)d_j=\\ & \underset{x_a}{\overset{x_b}{\int }}{\varphi }_{i}qdx+{\varphi }_i(x_a)Q_2+{\varphi }_i(x_b)Q_5+\frac{d{\varphi }_i}{dx}(x_a)Q_3+\frac{d{\varphi }_i}{dx}(x_b)Q_6\end{array}$

We can write the above as

$\sum _{j=1}^2{K}_{ij}^{21}{u}_j+\sum _{j=1}^4{K}_{ij}^{22}{d}_j=F_i^2$

where $i=1,2,3,4$ and

$\begin{array}{cc}{K}_{ij}^{21}& =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\varphi }_i}{dx}\frac{d{\psi }_j}{dx}dx\\ K_{ij}^{22}& =\underset{x_a}{\overset{x_b}{\int }}\{\frac{1}{2}{A}_{xx}[\frac{d{u}_0}{dx}+{\left(\frac{d{w}_0}{dx}\right)}^2]\frac{d{\varphi }_i}{dx}\frac{d{\varphi }_j}{dx}+D_{xx}\frac{d^2{\varphi }_i}{d{x}^2}\frac{d^2{\varphi }_j}{d{x}^2}\}dx\\ F_i^2& =\underset{x_a}{\overset{x_b}{\int }}{\varphi }_{i}qdx+{\varphi }_i(x_a)Q_2+{\varphi }_i(x_b)Q_5+\frac{d{\varphi }_i}{dx}(x_a)Q_3+\frac{d{\varphi }_i}{dx}(x_b)Q_6.\end{array}$

This gives us a symmetric stiffness matrix.

Part 6

Derive the element tangent stiffness matrix for the element.

Equation (21) can be written as

$𝐊(𝐔)𝐔=𝐅$

where

$\begin{array}{cc}{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{array}$

The residual is

$𝐑=𝐊𝐔-𝐅.$

For Newton iterations, we use the algorithm

${𝐔}^{r+1}={𝐔}^r-({𝐓}^r{)}^{-1}{𝐑}^r$

where the tangent stiffness matrix is given by

${𝐓}^r=\frac{\partial {𝐑}^r}{\partial 𝐔}.$

The coefficients of the tangent stiffness matrix are given by

$T_{ij}=\frac{\partial R_i}{\partial U_j},i=1\dots 6,j=1\dots 6.$

Recall that the finite element system of equations can be written as

$\begin{array}{cc}\sum _{p=1}^2{K}_{mp}^{11}{u}_p& +\sum _{q=1}^4{K}_{mq}^{12}{d}_q=F_m^1,m=1,2\\ \sum _{p=1}^2{K}_{np}^{21}{u}_p& +\sum _{q=1}^4{K}_{nq}^{22}{d}_q=F_n^2,n=1,2,3,4\end{array}$

where the subscripts have been changed to avoid confusion.

Therefore, the residuals are

$\begin{array}{cc}{R}_m^1& =\sum _{p=1}^2{K}_{mp}^{11}{u}_p+\sum _{q=1}^4{K}_{mq}^{12}{d}_q-F_m^1,m=1,2\\ R_n^2& =\sum _{p=1}^2{K}_{np}^{21}{u}_p+\sum _{q=1}^4{K}_{nq}^{22}{d}_q-F_n^2,n=1,2,3,4.\end{array}$

The derivatives of the residuals with respect to the generalized displacements are

$\begin{array}{cccc}{T}_{mk}^{11}=\frac{\partial R_m^1}{\partial u_k}& =\sum _{p=1}^2\frac{\partial }{\partial u_k}(K_{mp}^{11}{u}_p)+\sum _{q=1}^4\frac{\partial }{\partial u_k}(K_{mq}^{12}{d}_q),& & m=1,2;k=1,2\\ T_{ml}^{12}=\frac{\partial R_m^1}{\partial d_l}& =\sum _{p=1}^2\frac{\partial }{\partial d_l}(K_{mp}^{11}{u}_p)+\sum _{q=1}^4\frac{\partial }{\partial d_l}(K_{mq}^{12}{d}_q),& & m=1,2;l=1,2,3,4\\ T_{nk}^{21}=\frac{\partial R_n^2}{\partial u_k}& =\sum _{p=1}^2\frac{\partial }{\partial u_k}(K_{np}^{21}{u}_p)+\sum _{q=1}^4\frac{\partial }{\partial u_k}(K_{nq}^{22}{d}_q),& & n=1,2,3,4;k=1,2\\ T_{nl}^{22}=\frac{\partial R_n^2}{\partial d_l}& =\sum _{p=1}^2\frac{\partial }{\partial d_l}(K_{np}^{21}{u}_p)+\sum _{q=1}^4\frac{\partial }{\partial d_l}(K_{nq}^{22}{d}_q),& & n=1,2,3,4;l=1,2,3,4.\end{array}$

Differentiating, we get

$\begin{array}{cccc}{T}_{mk}^{11}& =\sum _{p=1}^2(u_p\frac{\partial K_{mp}^{11}}{\partial u_k}+K_{mp}^{11}\frac{\partial u_p}{\partial u_k})+\sum _{q=1}^4(d_q\frac{\partial K_{mq}^{12}}{\partial u_k}+K_{mq}^{12}\frac{\partial d_q}{\partial u_k}),& & m=1,2;k=1,2\\ T_{ml}^{12}& =\sum _{p=1}^2(u_p\frac{\partial K_{mp}^{11}}{\partial d_l}+K_{mp}^{11}\frac{\partial u_p}{\partial d_l})+\sum _{q=1}^4(d_q\frac{\partial K_{mq}^{12}}{\partial d_l}+K_{mq}^{12}\frac{\partial d_q}{\partial d_l}),& & m=1,2;l=1,2,3,4\\ T_{nk}^{21}& =\sum _{p=1}^2(u_p\frac{\partial K_{np}^{21}}{\partial u_k}+K_{np}^{21}\frac{\partial u_p}{\partial u_k})+\sum _{q=1}^4(d_q\frac{\partial K_{nq}^{22}}{\partial u_k}+K_{nq}^{22}\frac{\partial d_q}{\partial u_k}),& & n=1,2,3,4;k=1,2\\ T_{nl}^{22}& =\sum _{p=1}^2(u_p\frac{\partial K_{np}^{21}}{\partial d_l}+K_{np}^{21}\frac{\partial u_p}{\partial d_l})+\sum _{q=1}^4(d_q\frac{\partial K_{nq}^{22}}{\partial d_l}+K_{nq}^{22}\frac{\partial d_q}{\partial d_l}),& & n=1,2,3,4;l=1,2,3,4.\end{array}$

These equations can therefore be written as

$\begin{array}{cccc}{T}_{mk}^{11}& =K_{mk}^{11}+\sum _{p=1}^2{u}_p\frac{\partial K_{mp}^{11}}{\partial u_k}+\sum _{q=1}^4{d}_q\frac{\partial K_{mq}^{12}}{\partial u_k},& & m=1,2;k=1,2\\ T_{ml}^{12}& =K_{ml}^{12}+\sum _{p=1}^2{u}_p\frac{\partial K_{mp}^{11}}{\partial d_l}+\sum _{q=1}^4{d}_q\frac{\partial K_{mq}^{12}}{\partial d_l},& & m=1,2;l=1,2,3,4\\ T_{nk}^{21}& =K_{nk}^{21}+\sum _{p=1}^2{u}_p\frac{\partial K_{np}^{21}}{\partial u_k}+\sum _{q=1}^4{d}_q\frac{\partial K_{nq}^{22}}{\partial u_k},& & n=1,2,3,4;k=1,2\\ T_{nl}^{22}& =K_{nl}^{22}+\sum _{p=1}^2{u}_p\frac{\partial K_{np}^{21}}{\partial d_l}+\sum _{q=1}^4{d}_q\frac{\partial K_{nq}^{22}}{\partial d_l},& & n=1,2,3,4;l=1,2,3,4.\end{array}$

Now, the coefficients of ${𝐊}^{11}$, ${𝐊}^{12}$, and ${𝐊}^{21}$ of the symmetric stiffness matrix are independent of $u_1$ and $u_2$. Also, the terms of ${𝐊}^{11}$ are independent of the all the generalized displacements. Therefore, the above equations reduce to

$\begin{array}{cccc}{T}_{mk}^{11}& =K_{mk}^{11},& & m=1,2;k=1,2\\ T_{ml}^{12}& =K_{ml}^{12}+\sum _{q=1}^4{d}_q\frac{\partial K_{mq}^{12}}{\partial d_l},& & m=1,2;l=1,2,3,4\\ T_{nk}^{21}& =K_{nk}^{21}+\sum _{q=1}^4{d}_q\frac{\partial K_{nq}^{22}}{\partial u_k},& & n=1,2,3,4;k=1,2\\ T_{nl}^{22}& =K_{nl}^{22}+\sum _{p=1}^2{u}_p\frac{\partial K_{np}^{21}}{\partial d_l}+\sum _{q=1}^4{d}_q\frac{\partial K_{nq}^{22}}{\partial d_l},& & n=1,2,3,4;l=1,2,3,4.\end{array}$

Consider the coefficients of ${𝐓}^{12}$:

$T_{ml}^{12}=K_{ml}^{12}+\sum _{q=1}^4{d}_q\frac{\partial K_{mq}^{12}}{\partial d_l},m=1,2;l=1,2,3,4.$

From our previous derivation, we have

$K_{mq}^{12}=\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\psi }_m}{dx}\frac{d{\varphi }_q}{dx}dx.$

Therefore,

$\begin{array}{cc}\frac{\partial K_{mq}^{12}}{\partial d_l}& =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left[A_{xx}\frac{\partial }{\partial d_l}\left(\frac{d{w}_0}{dx}\right)\right]\frac{d{\psi }_m}{dx}\frac{d{\varphi }_q}{dx}dx\\ & =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left[A_{xx}\left(\sum _{T=1}^4\frac{\partial d_T}{\partial d_l}\frac{d{\varphi }_T}{dx}\right)\right]\frac{d{\psi }_m}{dx}\frac{d{\varphi }_q}{dx}dx\\ & =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left[A_{xx}\frac{d{\varphi }_l}{dx}\right]\frac{d{\psi }_m}{dx}\frac{d{\varphi }_q}{dx}dx\end{array}$

The tangent stiffness matrix coefficients are therefore

$\begin{array}{cc}{T}_{ml}^{12}& =K_{ml}^{12}+\sum _{q=1}^4{d}_q\left\{\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left[A_{xx}\frac{d{\varphi }_l}{dx}\right]\frac{d{\psi }_m}{dx}\frac{d{\varphi }_q}{dx}dx\right\}m=1,2;l=1,2,3,4.\\ & =K_{ml}^{12}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\varphi }_l}{dx}\frac{d{\psi }_m}{dx}\left(\sum _{q=1}^4{d}_q\frac{d{\varphi }_q}{dx}\right)dx\\ & =K_{ml}^{12}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\varphi }_l}{dx}\frac{d{\psi }_m}{dx}\frac{d{w}_0}{dx}dx\\ & =2K_{ml}^{12}.\end{array}$

Next, {consider the coefficients of ${𝐓}^{21}$:

$T_{nk}^{21}=K_{nk}^{21}+\sum _{q=1}^4{d}_q\frac{\partial K_{nq}^{22}}{\partial u_k}.$

The coefficients of ${𝐊}^{22}$ are

$K_{nq}^{22}=\underset{x_a}{\overset{x_b}{\int }}\{\frac{1}{2}{A}_{xx}[\frac{d{u}_0}{dx}+{\left(\frac{d{w}_0}{dx}\right)}^2]\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}+D_{xx}\frac{d^2{\varphi }_n}{d{x}^2}\frac{d^2{\varphi }_q}{d{x}^2}\}dx.$

Therefore, the derivatives are

$\begin{array}{cc}\frac{\partial K_{nq}^{22}}{\partial u_k}& =\underset{x_a}{\overset{x_b}{\int }}\frac{1}{2}{A}_{xx}[\sum _{T=1}^2\frac{\partial u_T}{\partial u_k}\frac{d{\psi }_T}{dx}+2\frac{d{w}_0}{dx}\left(\sum _{T=1}^4\frac{\partial d_T}{\partial u_k}\frac{d{\varphi }_T}{dx}\right)]\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}dx\\ & =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\psi }_k}{dx}\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}dx\\ \end{array}$

Therefore the coefficients of ${𝐓}^{21}$ are

$\begin{array}{cc}{T}_{nk}^{21}& =K_{nk}^{21}+\sum _{q=1}^4{d}_q\left(\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\psi }_k}{dx}\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}dx\right)\\ & =K_{nk}^{21}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\psi }_k}{dx}\frac{d{\varphi }_n}{dx}\left(\sum _{q=1}^4{d}_q\frac{d{\varphi }_q}{dx}\right)dx\\ & =K_{nk}^{21}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\psi }_k}{dx}\frac{d{\varphi }_n}{dx}\frac{d{w}_0}{dx}dx\\ & =2K_{nk}^{21}.\end{array}$

Finally, for the ${𝐓}^{22}$ coefficients, we start with

$T_{nl}^{22}=K_{nl}^{22}+\sum _{p=1}^2{u}_p\frac{\partial K_{np}^{21}}{\partial d_l}+\sum _{q=1}^4{d}_q\frac{\partial K_{nq}^{22}}{\partial d_l},n=1,2,3,4;l=1,2,3,4$

and plug in the derivatives of the stiffness matrix coefficients

$\begin{array}{cc}{K}_{np}^{21}& =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left(A_{xx}\frac{d{w}_0}{dx}\right)\frac{d{\varphi }_n}{dx}\frac{d{\psi }_p}{dx}dx\\ K_{nq}^{22}& =\underset{x_a}{\overset{x_b}{\int }}\{\frac{1}{2}{A}_{xx}[\frac{d{u}_0}{dx}+{\left(\frac{d{w}_0}{dx}\right)}^2]\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}+D_{xx}\frac{d^2{\varphi }_n}{d{x}^2}\frac{d^2{\varphi }_q}{d{x}^2}\}dx.\end{array}$ \\

The derivatives are

$\begin{array}{cc}\frac{\partial K_{np}^{21}}{\partial d_l}& =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}\left[A_{xx}\frac{\partial }{\partial d_l}\left(\frac{d{w}_0}{dx}\right)\right]\frac{d{\varphi }_n}{dx}\frac{d{\psi }_p}{dx}dx\\ & =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\left[\sum _{T=1}^4\frac{\partial d_T}{\partial d_l}\frac{d{\varphi }_T}{dx}\right]\frac{d{\varphi }_n}{dx}\frac{d{\psi }_p}{dx}dx\\ & =\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\frac{d{\psi }_p}{dx}dx\end{array}$

and

$\begin{array}{cc}\frac{\partial K_{nq}^{22}}{\partial d_l}& =\underset{x_a}{\overset{x_b}{\int }}\frac{1}{2}{A}_{xx}[\sum _{T=1}^2\frac{\partial u_T}{\partial d_l}\frac{d{\psi }_T}{dx}+2\frac{d{w}_0}{dx}\left(\sum _{T=1}^4\frac{\partial d_T}{\partial d_l}\frac{d{\varphi }_T}{dx}\right)]\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}dx\\ & =\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{w}_0}{dx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}dx\end{array}$

Therefore, the coefficients of the tangent stiffness matrix can be written as

$\begin{array}{cc}{T}_{nl}^{22}& =K_{nl}^{22}+\sum _{p=1}^2{u}_p\left[\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\frac{d{\psi }_p}{dx}dx\right]+\sum _{q=1}^4{d}_q\left[\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{w}_0}{dx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\frac{d{\varphi }_q}{dx}dx\right]\\ & =K_{nl}^{22}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\left(\sum _{p=1}^2{u}_p\frac{d{\psi }_p}{dx}\right)dx+\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{w}_0}{dx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\left(\sum _{q=1}^4{d}_q\frac{d{\varphi }_q}{dx}\right)dx\\ & =K_{nl}^{22}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\frac{d{u}_0}{dx}dx+\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}\frac{d{w}_0}{dx}\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}\frac{d{w}_0}{dx}dx\\ & =K_{nl}^{22}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}[\frac{d{u}_0}{dx}+2{\left(\frac{d{w}_0}{dx}\right)}^2]\frac{d{\varphi }_l}{dx}\frac{d{\varphi }_n}{dx}dx.\end{array}$

The final expressions for the tangent stiffness matrix terms are

$\begin{array}{cccc}{T}_{ij}^{11}& =K_{ij}^{11},& & i=1,2;j=1,2\\ T_{ij}^{12}& =2K_{ij}^{12},& & i=1,2;j=1,2,3,4\\ T_{ij}^{21}& =2K_{ij}^{21},& & i=1,2,3,4;j=1,2\\ T_{ij}^{22}& =K_{ij}^{22}+\frac{1}{2}\underset{x_a}{\overset{x_b}{\int }}{A}_{xx}[\frac{d{u}_0}{dx}+2{\left(\frac{d{w}_0}{dx}\right)}^2]\frac{d{\varphi }_i}{dx}\frac{d{\varphi }_j}{dx}dx& & i=1,2,3,4;j=1,2,3,4.\end{array}$