Modeling the St. Louis Arch

Posted on 2019-10-28 in how-to • 4 min read

Learning By Doing

Last year over Spring Break, we visited St. Louis. Our trip included a stop at the arch, and we managed to take a successful trip to the observation deck at the top. Shortly thereafter, I happened across some information on the National Park Service page regarding the geometry of the arch. Being a civil engineer with an acute interest in geometry, I couldn’t resist tinkering around with these equations as a means of diving into computational geometry.

Research

A good first step before diving into any scientific challenge is to perform research. In doing so, I found a number of resources in addition to the official NPS website:

Gateway to Mathematics Equations of the St. Louis Arch

How the Gateway Arch Got its Shape

Owner’s Manual For the Gateway Arch

This initial research led to two important initial findings:

  1. The geometry of the arch can indeed be specified mathematically.
  2. The primary form of the arch is derived from the shape of a flattened catenary curve.

Catenary curves

Wikipedia is a good starting point for exploring the computational geometry needed to accurately model the arch. The basic equation for a catenary involves a hyperbolic cosine (cosh) and a single constant that defines the aspect ratio, or flattening factor, of the resulting curve.

In [7]:
import math

import matplotlib.pyplot as plt
plt.rcParams['figure.figsize'] = [11, 8.5]
import numpy as np
import mpmath 

def catenary(x, a):
    y = a * np.cosh(x / a)
    return y

x = np.arange(-4, 4, 0.1)
y1 = catenary(x, 1.0)
y2 = catenary(x, 2.0) 
yhalf = catenary(x, 0.5)

plt.plot(x, y1, label='a = 1.0')
plt.plot(x, y2, label='a = 2.0')
plt.plot(x, yhalf, label='a = 0.5')
axes = plt.gca()
axes.set_ylim(0, 6)
plt.title('Fig. 1: Sample Catenary Curves')
plt.legend()
plt.show()

Flattened Catenary

As mentioned in the reference documents, the shape of the arch is a modified, or flattened, catenary curve. This means that there are two constants needed in the catenary equation to properly scale and proportion the initial curve that will form the basis of the arch geometry.

Constants from NPS web site

To develop a mathematical model for the arch geometry, the limits of X and Y need to be developed. Additionally, the arch is a 3-dimensional shape created by a series of equilateral triangles that are larger at the bottom and smaller at the top. A flattened catenary equation provides the location of the centroid of each triangle. Two additional values are required in order to set the “width” or “thickness” of the arch at a given location.

These four values (width at base, height at apex, maximum thickness, and minimum thickness) are listed on the NPS website as follows:

In [8]:
fc = 625.0925 # Maximium height of centroid
Qb = 1262.6651 # Maximum cross-sectional area at arch base
Qt = 125.1406 # Minimum cross-sectional area at arch top
L = 299.2239 # Half of centroid at arch base

Equations from NPS web site

As noted previously, the arch is a flattened catenary and therefore requires two constants. These constants were influenced by the minimum and maximum cross-sectional areas as noted above. The constants are listed on the NPS website as follows:

In [9]:
A = fc / ((Qb / Qt) - 1)
C = math.acosh(Qb / Qt)

print(f"A = {A:.4f}")
print(f"C = {C:.4f}")

A = 68.7673
C = 3.0022

Interestingly, there is a small rounding difference in the fourth decimal place of the calculated value for A when comparing it to the published value of 68.7672.

Next, we develop a function to calculate the value of X at any provided Y, where Y is the vertical axis. However, there is an important adjustment required: A typical catenary curve “opens up” as Y increases upwards as shown in Fig. 1.

This would be a very unstable arrangement for the arch. Therefore, the equations provided in the construction blueprints inverted the direction of the Y axis so that the catenary would have maximum width at ground level. In other words, the positive direction of Y is towards the ground.

In construction, we typically think in terms of elevation as a distance above ground level. Therefore, we will define a “helper function” to convert Y values to elevation before plugging that value into the function for X.

Note that because we are solving for X at a given Y value instead of the other way around as shown earlier, the arch equation uses the inverse hyperbolic cosine.

In [10]:
def El(Y):
    return fc - Y

def X(Y):
    return (L / C) * (math.acosh(1 + El(Y)/A))

Let’s test our function. We know that the maximum height of the arch is given by fc. At Y = fc, X should be zero. Correspondingly, at Y = 0, X should be at a maximum, or roughly equal to L.

In [11]:
print(f"Value of X at Y = fc is {X(fc):.4f}")
print(f"Value of X at Y = 0 is {X(0):.4f}")
print(f"L = {L}")

Value of X at Y = fc is 0.0000
Value of X at Y = 0 is 299.2239
L = 299.2239

Let’s plot our function and see how things are looking. We will need to vectorize the function to apply it over an array of elevation values. We also need to plot a second, mirrored version so that we can get the whole picture.

In [12]:
elevations = np.arange(0, fc, 0.1)

X_vectorized = np.vectorize(X, otypes=[np.float])

centroid = X_vectorized(elevations) 

plt.plot(centroid, elevations)
plt.plot(-centroid, elevations)
axes = plt.gca()
axes.set_aspect('equal')
plt.title('Fig. 2: Centroid of Gateway Arch')
plt.show()

Interestingly, this plot mirrors the construction sequence where the two legs of the arch were constructed individually and then met at the top to form a single structure.

Conclusion

We’ve gotten a good start on developing our model via computational geometry. However, we certainly can’t call this a 4D VDC model yet. The next post will look at adding the triangular cross sections in order to calculate the inner (intrados) and outer (extrados) edges of the arch.

python