Intended learning outcomes

• Derive the strong form of the equilibrium equations for a solid

• Apply appropriate boundary conditions to the boundary value problem

The boundary value problem (BVP)

Our goal is to investigate a structure subjected to external loads. We might, for example, be interested in the resulting deformations or the internal stresses. One way to think about this, is that applying known displacements to the surface causes stresses, $\boldsymbol{ \sigma}$, inside the body. For equilibrium to hold, we require that these stresses fulfill

$$\begin{aligned} \boldsymbol{ \sigma}\cdot\underline{\boldsymbol{ \nabla}} + \underline{\boldsymbol{ b}} &= \underline{\boldsymbol{ 0}} \end{aligned}$$

inside the body subjected to the body load $\underline{\boldsymbol{ b}}$ [Force/Volume], and that the traction $\underline{\boldsymbol{ t}}$ is

$$\begin{aligned} \underline{\boldsymbol{ t}} &= \underline{\boldsymbol{ n}}\cdot\boldsymbol{ \sigma} \end{aligned}$$

on the boundary, where $\underline{\boldsymbol{ n}}$ is the normal vector.

How high the stresses become for given displacements, depends on the material behavior. That is, how $\boldsymbol{ \sigma}$ depends on the strain $\boldsymbol{ \epsilon}$. For linear elasticity, we have that

$$\begin{aligned} \boldsymbol{ \sigma}(\boldsymbol{ \epsilon}) = \textbf{\textsf{ E}}:\boldsymbol{ \epsilon}, \quad \text{(Linear elasticity)} \end{aligned}$$

We thus need to determine the strains, $\boldsymbol{ \epsilon}$, inside the body to calculate the stresses. And the strains depend on the displacements $\underline{\boldsymbol{ u}}$ as

$$\begin{aligned} \boldsymbol{ \epsilon}(\underline{\boldsymbol{ u}}) = \left[ \underline{\boldsymbol{ u}}\otimes\underline{\boldsymbol{ \nabla}}\right]^{ \mathrm{sym} } \end{aligned}$$

We now have the differential equations required to solve the boundary value problem, namely

$$\begin{aligned} \boldsymbol{ \sigma}(\boldsymbol{ \epsilon}(\underline{\boldsymbol{ u}}(\underline{\boldsymbol{ x}})))\cdot\underline{\boldsymbol{ \nabla}} + \underline{\boldsymbol{ b}} &= \underline{\boldsymbol{ 0}},\phantom{(\underline{\boldsymbol{ x}})} \quad \underline{\boldsymbol{ x}} \in \Omega \\ \underline{\boldsymbol{ u}}(x) &= \underline{\boldsymbol{ u}}_{ \mathrm{ D} }(\underline{\boldsymbol{ x}}), \quad \underline{\boldsymbol{ x}} \in \Gamma_{ \mathrm{ D} } \\ \underline{\boldsymbol{ n}} \cdot \boldsymbol{ \sigma}(\boldsymbol{ \epsilon}(\underline{\boldsymbol{ u}}(\underline{\boldsymbol{ x}}))) &= \underline{\boldsymbol{ t}}_{ \mathrm{ N} }(\underline{\boldsymbol{ x}}), \quad \underline{\boldsymbol{ x}} \in \Gamma_{ \mathrm{ N} } \end{aligned}$$

where the first equation defines the differential equation inside the body (light orange). The second equation give the so-called Dirichlet boundary conditions (black, $\Gamma_{ \mathrm{ d} }$), here we already know the value of the displacements. And finally, the part of the boundary with unknown displacements is called the Neuman (also natural) boundary (blue, $\Gamma_{ \mathrm{ N} }$). Here, we know the traction. Often, a majority of the boundary has zero load but we know it is zero!

The problem is - it is not possible to determine the function $\underline{\boldsymbol{ u}}(\underline{\boldsymbol{ x}})$ for most cases. While it is possible in a few special cases, what to do for the remaining problems? Often, we use the finite element method. But now let's first start with a case where we can find the solution analytically.

Analytical solution to the BVP (uniaxial stress)

Consider the case of cylinder subjected to uniaxial tension, with an isotropic material:

On the blue edges we have Neumann (natural) boundary conditions. We don't know the displacements, but we know the traction. On the right side ($x_1=L$), this traction is $[F/A]\underline{\boldsymbol{ e}}_{ 1}$, where $A$ is the cross-sectional area and on the sides of the cylinder the traction is zero. On the left side, we actually have mixed boundary conditions. We know that the traction components in the $x_2$-$x_3$ plane are zero (Neumann), and that the displacement is zero in the $\underline{\boldsymbol{ e}}_{ 1}$ direction (Dirichlet).

To solve the differential equation, we make the following ansatz $$\begin{aligned} \underline{\boldsymbol{ u}}(\underline{\boldsymbol{ x}}) = a_1 x_1 \underline{\boldsymbol{ e}}_{ 1} + a_2 x_2 \underline{\boldsymbol{ e}}_{ 2} + a_3 x_3 \underline{\boldsymbol{ e}}_{ 3} \end{aligned}$$

From this ansatz, we obtain the constant strains as

$$\begin{aligned} \boldsymbol{ \epsilon} = a_1 \underline{\boldsymbol{ e}}_{ 1}\otimes\underline{\boldsymbol{ e}}_{ 1} + a_2 \underline{\boldsymbol{ e}}_{ 2}\otimes\underline{\boldsymbol{ e}}_{ 2} + a_3 \underline{\boldsymbol{ e}}_{ 3}\otimes\underline{\boldsymbol{ e}}_{ 3} \end{aligned}$$

For linear isotropic elasticity, we have $\boldsymbol{ \sigma} = 2G \boldsymbol{ \epsilon}^{ \mathrm{dev} } + 3K\boldsymbol{ \epsilon}^{ \mathrm{sph} }$, with

$$\begin{aligned} \boldsymbol{ \epsilon}^{ \mathrm{sph} } &= \frac{a_1+a_2+a_3}{3} \boldsymbol{ I} \\ \boldsymbol{ \epsilon}^{ \mathrm{dev} } &= \frac{2a_1 - a_2 - a_3}{3}\underline{\boldsymbol{ e}}_{ 1}\otimes\underline{\boldsymbol{ e}}_{ 1} + \frac{2a_2 - a_1 - a_3}{3}\underline{\boldsymbol{ e}}_{ 2}\otimes\underline{\boldsymbol{ e}}_{ 2} + \frac{2a_3 - a_1 - a_2}{3}\underline{\boldsymbol{ e}}_{ 3}\otimes\underline{\boldsymbol{ e}}_{ 3} \end{aligned}$$

Checking the equilibrium by taking the divergence of $\boldsymbol{ \sigma}$, we see that this is zero as $\boldsymbol{ \sigma}$ is constant. Hence, the proposed ansatz fulfills the equilbrium equation. However, we must also fulfill the boundary conditions. At $x_1=L$, we have that $\underline{\boldsymbol{ t}}=[F/A]\underline{\boldsymbol{ e}}_{ 1}$, that is

$$\begin{aligned} \underline{\boldsymbol{ t}} &= \underline{\boldsymbol{ e}}_{ 1} \cdot \boldsymbol{ \sigma} = \underline{\boldsymbol{ e}}_{ 1} \cdot [2G \boldsymbol{ \epsilon}^{ \mathrm{dev} } + 3K\boldsymbol{ \epsilon}^{ \mathrm{sph} }] \\ &= 2G \frac{2a_1 - a_2 - a_3}{3} \underline{\boldsymbol{ e}}_{ 1} + K [a_1+a_2+a_3] \underline{\boldsymbol{ e}}_{ 1} = \frac{F}{A}\underline{\boldsymbol{ e}}_{ 1} \\ \frac{F}{A} &= 2G \frac{2a_1 - a_2 - a_3}{3} + K [a_1+a_2+a_3] \end{aligned}$$

Considering the boundary conditions on the side of the cylinder, e.g. at $x_2=r,\,x_3=0$ and $x_2=0,\,x_3=r$, we have that the traction is zero, e.g.

$$\begin{aligned} \underline{\boldsymbol{ t}} &= \underline{\boldsymbol{ e}}_{ 2} \cdot \boldsymbol{ \sigma} = 2G \frac{2a_2 - a_1 - a_3}{3} \underline{\boldsymbol{ e}}_{ 2} + K [a_1+a_2+a_3] \underline{\boldsymbol{ e}}_{ 2} = \underline{\boldsymbol{ 0}} \\ 0 &= 2G \frac{2a_2 - a_1 - a_3}{3} + K [a_1+a_2+a_3] \\ 0 &= 2G \frac{2a_3 - a_1 - a_2}{3} + K [a_1+a_2+a_3], \quad \mathrm{obtained\,the\,same\,way\,as\,above} \end{aligned}$$

And hence, we see from the two last equivalent equations where we just swapped $a_2$ and $a_3$ that $a_2=a_3$. Resulting in

$$\begin{aligned} \frac{F}{A} &= 4G \frac{a_1 - a_2}{3} + K [a_1+2a_2] \\ 0 &= 2G \frac{a_2 - a_1}{3} + K [a_1+2a_2] \\ a_2 &= -\frac{3K - 2G}{2G + 6K}a_1 \end{aligned}$$

To simplify, further calculation, we denote $\nu = [3K-2G]/[2G+6K]$ (This is the actual Poissons' ratio, as we can see from how it governs the relation between $a_1$ and $a_2$ along with Equation (7)). We then have

$$\begin{aligned} a_2 &= -\nu a_1 \\ \frac{F}{A} &= 4G \frac{a_1 + \nu a_1}{3} + K [a_1-2\nu a_1] \\ a_1 &= \frac{F}{A[4G[1+\nu] + K[1-2\nu]]} \end{aligned}$$

This expression is now simplified by denoting $E=4G[1+\nu] + K[1-2\nu]$ (This is the actual Young's modulus, as we will se later). And we have that

$$\begin{aligned} a_1 &= \frac{F}{AE} \\ a_2 &= -\nu \frac{F}{AE} = a_3 \end{aligned}$$

Yielding the following displacements and strains.

$$\begin{aligned} \underline{\boldsymbol{ u}}(\underline{\boldsymbol{ x}}) &= \frac{F}{AE} \left[ x_1 \underline{\boldsymbol{ e}}_{ 1} - \nu x_2 \underline{\boldsymbol{ e}}_{ 2} - \nu x_3 \underline{\boldsymbol{ e}}_{ 3}\right] \\ \boldsymbol{ \epsilon} &= \frac{F}{AE} \left[ \underline{\boldsymbol{ e}}_{ 1}\otimes\underline{\boldsymbol{ e}}_{ 1} - \nu [\underline{\boldsymbol{ e}}_{ 2}\otimes\underline{\boldsymbol{ e}}_{ 2} + \underline{\boldsymbol{ e}}_{ 3}\otimes\underline{\boldsymbol{ e}}_{ 3}]\right] \end{aligned}$$

Finally, if we insert the constitutive relationship $\boldsymbol{ \sigma}=2G \boldsymbol{ \epsilon}^{ \mathrm{dev} } + K\boldsymbol{ \epsilon}^{ \mathrm{sph} }$, and using our definitions of $E$ and $\nu$ from above, we obtain

$$\begin{aligned} \boldsymbol{ \sigma} &= \frac{F}{A} \underline{\boldsymbol{ e}}_{ 1}\otimes\underline{\boldsymbol{ e}}_{ 1} \end{aligned}$$

I.e. we see that the stress is indeed uniaxial. Furthermore, seeing combining Equation (14) with Equation (15), we see that $\sigma_{11}=E \epsilon_{11}$, showing that the $E$ we defined earlier is the actual Young's modulus.

In summary, for even this very simple case, determining the solution to the BVP in Equation (5) was quite a bit of effort. And as previously noted, in many cases no unique solution exists even if we include simplifications such as plane stress or strain.