In my ongoing quest to reinvent all the equations in all the engineering codes, I found myself idly wondering why the equations of shell thickness in API 650 are the way they are. You would think that the governing stress in a storage tank is the hoop stress, and the max hoop stress is at the bottom of the shell course where the hydrostatic pressure is the greatest. However if you look at the code you find the minimum thickness for a tank with diameter $D$ containing a fluid with specific gravity $SG$ and fluid height $H$:

$$ t = { { 4.9 \cdot SG \cdot D \cdot (H - x)} \over S } $$

Where $x$ is a distance somewhere up from the bottom of the shell course, to be determined either empirically as "about a foot" or through the complicated set of madness that is the variable design point method (an issue for another day).

### a first approach

At first glance you may be wondering "what is this crazy equation?" Lets back up and start with the equation for hoop stress:

$$ \sigma = { { p \cdot D} \over 2t} $$

The hydrostatic pressure at some point $x$, measured from the bottom of the shell course is:

$$ p = \rho g (H-x) $$

So the hoop stress at that point is:

$$ \sigma = {{ \rho \cdot g \cdot D \cdot (H-x) } \over 2t} = {{4.9 \cdot SG \cdot D \cdot (H-x)} \over t} $$

$$ t = {{ 4.9 \cdot SG \cdot D \cdot (H-x)} \over S} $$

Where I've substituted in the specific gravity and the unit conversions such that the stress is in MPa and the thickness in mm instead of the SI base units of Pa and m.

So that half explains that, API 650 is calculating the minimum thickness needed to withstand the hoop stress at a height $x$. But why is it some $x$ and not the bottom of the shell course, where the hoop stress would naively be higher?

### displacement in a thin shell

There are a couple of different ways of addressing this question:

- Pull out Kirchoff-Love plate theory and solve some ODEs
- Pull out Roark's and look at a superposition of stresses

I am going to be lazy and look at Roark's, i.e. *Roark's Formulas for Stress and Strain, 7th Edition*. The nice thing about Roark's is that it collates solutions to the various stress equations for common loadings, so basically someone else has solved those differential equations for you. What remains is to perform superposition to get the particular set of loads that your system experiences.

I'm going to look at a hypothetical shell course with hydrostatic pressure, supported at the bottom in such a way that the radial deformation is zero and the rotation is zero. The classic situation would be if it was welded to the floor plates.

This is a superposition of 3 situations from Roark's:

- Situation 1(d) from Table 13.1
- Situation 8 from Table 13.2
- Situation 10 from Table 13.2

Note$D$ is not the tank diameter from here on out, it is the bending rigidity:$$ D = \frac{Et^3}{12 (1 - \nu^2)} $$

I am just assuming that since storage tanks are typically large relative to the thickness of the shell, that $\lambda H > 6$, which allows me to use 8 and 10 from Table 13.2

The plan of attack is this:

- Find an expression for the radial deformation at the bottom
- Find an expression for the rotation at the bottom
- Solve (1) and (2) for the load $V_{o}$ and $M_{o}$
- Use the values from (3) to get an expression for the radial deformation as a function of $y$

First off the radial deformation at the bottom is zero, by assumption, and it is the sum of the contributions from each situation:

Situation 1(d) from Table 13.1

$$ \Delta R = \frac{q_{0} R^2}{Et} $$

Situation 8 from Table 13.2

$$ \Delta R = - \frac{V_{o}}{2D \lambda^3} $$

Situation 10 from Table 13.2

$$ \Delta R = \frac{M_{o}}{2D \lambda^2} $$

With $\lambda$ given by:

$$ \lambda = \left[ \frac{3 (1 - \nu^2)}{(Rt)^2} \right]^{\frac{1}{4}} $$

Anywho, the final expression for the radial deformation at the bottom is (dropping the subscript "o")

$$ \Delta R = \frac{q_{0} R^2}{Et} - \frac{V}{2D \lambda^3} + \frac{M}{2D \lambda^2} = 0 $$

Similarly the rotation about the bottom can be summed,

$$ \psi = \frac{V}{2D \lambda^2} - \frac{M}{D \lambda} = 0 $$

This is a system of linear equations with 2 equations and 2 unknowns, $V$ and $M$, solving $\psi$ for $M$

$$ M = \frac{V}{2 \lambda} $$

Subbing into the equation for $\Delta R$:

$$ V = 4D \lambda^3 \frac{q_{0} R^2}{Et} $$

and

$$ M = 2D \lambda^2 \frac{q_{0} R^2}{Et} $$

The radial deflection at any point $x$ along the shell course is the superposition of the radial deflection due to the hydrostatic pressure, i.e. $q_{0} = \rho g H$ and the reactions at the bottom:

$$ \Delta R (x) = \frac{\rho g R^2}{Et} \left[ H - x - H \exp{(-\lambda x)} \left( \cos{\lambda x} + \sin{\lambda x} \right) \right] $$

The hoop stress is, by Hooke's law, proportional to the relative radial deflection:

$$ \sigma (x) = E \frac{\Delta R}{R} = \frac{\rho g R}{t} \left[ H - x - H \exp{(-\lambda x)} \left( \cos{\lambda x} + \sin{\lambda x} \right) \right] $$

### considering the result

Just glaring at the equation tells us something interesting. First off the hoop stress is actually quite small immediately adjacent to the shell bottom, due to the reaction forces there, but the impact of those forces decays as we move up the shell course -- this is the impact of $\exp{ (-\lambda x) }$

We would expect that there is a threshold where the hoop stress due to hydrostatic pressure comes to dominate and the hoop stress would then be linear with distance up to the top of the tank. The max. hoop stress would then be approximately:

$$ \sigma = \frac{\rho g R}{t} \left( H - x \right) = \frac{\rho g D}{2t} \left( H - x \right) = {{4.9 \cdot SG \cdot D \cdot (H-x)} \over t}$$

Where $x$ is the height where $\exp{ (-\lambda x) }$ is *basically zero* (waves hands in the air)

### exploring the result in code

Below are expressions for the hoop stress and radial displacement implemented in python.

```
from math import sin, cos, exp
def lmbda(nu,R,t):
return (3*(1-nu**2)/(R*t)**2)**(0.25)
def hoop_stress(rho,H,R,t):
''' returns a function for calculating the hoop stress at a point x from the base'''
nu = 0.3
g = 9.81 # m/s^2
l = lmbda(nu,R,t)
assert (l*H > 6), 'load criteria'
return lambda x: ((rho*g*R)/t)*(H-x - H*exp(-l*x)*(cos(l*x)+sin(l*x)))
def radial_disp(rho,H,R,t):
''' returns a function for calculating the radial displacement from a point x from the base'''
nu = 0.3
E = 200*1e9 #GPa
g = 9.81 #m/s^2
l = lmbda(nu,R,t)
assert (l*H > 6), 'load criteria'
return lambda x: ((rho*g*R*R)/(E*t))*(H-x - H*exp(-l*x)*(cos(l*x)+sin(l*x)))
```

We can plot what the hoop stress is doing for an arbitrary tank with a height of 4m, diameter of 2m, and a shell thickness of 5mm (full of water with a density of 1000kg/m3)

We can easily find the max of this function using scipy

```
import scipy.optimize as opt
max_val = opt.minimize_scalar(lambda x: -sigma(x), bounds=[0,0.5], method='bounded')
# message: 'Solution found.'
# x: 0.16517897498280573
# status: 0
# fun: -7855891.2806365527
# nfev: 12
# success: True
```

Comparing the *real* answer to the naive approach of $\sigma = \frac{\rho g R}{t} \left( H - x \right)$

```
real_sigma = sigma(max_val.x)
naive_sigma = (1000*9.81*1.0*(4.0-max_val.x))/(0.005)
per_diff = 100 * abs(real_sigma - naive_sigma)/real_sigma
print("Real answer {:f} MPa".format(real_sigma*1e-6))
print("Naive answer {:f} MPa".format(naive_sigma*1e-6))
print("% Difference {:f} %".format(per_diff))
# Real answer 7.855891 MPa
# Naive answer 7.523919 MPa
# % Difference 4.225777 %
```

Well that's interesting. I would have expected the naive answer to be *higher* than the "real" answer. If we plot just the component of the hoop stress from the reaction forces we see that while it decays towards zero it dips below and continues to oscillate around zero. This isn't too shocking since there are those cos and sin functions in there.

Regardless of the specifics, the error is small in making the naive assumption and it is more workable for finding the required shell thickness since $\lambda$ depends upon thickness. When all the safety factors are included the slight error in not pinpointing the *exact* location of highest stress is negligible.

I assume that the Variable Design Point method does something similar but likely uses a more comprehensive model of the tank. For example I ignored the weight of the shell entirely as well as any impact the girth welds have on the stiffness of the shell. Once all this is included the hoop stress becomes a function of $x$ and $t$, which makes solving for $t$ all that more annoying.

As usual the ipython notebook is on github