# How Lighthouses Work 2, Exercise 3: why a tower is wider at the base
# An imagined tower 30 m high, 10 m across at the base and 6 m at the top, solid
# granite (2650 kg/m3 assumed). Find the average compressive stress at heights 0, 10
# and 20 m, compare its weight with a plain 10 m cylinder of the same height, and
# the toy dovetail resistance of one block against a ring of 10 (mu 0.5, block 3 t,
# 5 t of masonry resting on it).
import math

# An imagined tower (not a real one): height, base and top diameter, granite density
H, D_BASE, D_TOP = 30.0, 10.0, 6.0
RHO, G = 2650.0, 9.81

def diameter_at(y):
    return D_BASE + (D_TOP - D_BASE) * y / H

def area(d):
    return math.pi * d ** 2 / 4

def weight_above(y, steps=2000):
    dy = (H - y) / steps
    return sum(area(diameter_at(y + (i + 0.5) * dy)) * dy for i in range(steps)) * RHO * G

for y in (0, 10, 20):
    print(y, round(weight_above(y) / area(diameter_at(y)) / 1e3), "kPa")

# the same height as a plain cylinder of base diameter
cyl_weight = area(D_BASE) * H * RHO * G
print(round(weight_above(0) / cyl_weight, 2))
print(round(cyl_weight / area(D_BASE) / 1e3), "kPa (cylinder, at the base)")

# dovetail toy: ring of n blocks
MU, BLOCK_T, LOAD_T = 0.5, 3.0, 5.0
single = MU * (BLOCK_T + LOAD_T) * 1000 * G / 1000
print(round(single), "kN for one block; ring of 10:", round(single * 10), "kN")

# Output when run:
# 0 510 kPa
# 10 376 kPa
# 20 216 kPa
# 0.65
# 780 kPa (cylinder, at the base)
# 39 kN for one block; ring of 10: 392 kN
#
# Lines 1-3: stress rises from 216 kPa at 20 m up to 510 kPa at the base, because the
# base carries all the stone above it. The tapered form weighs 0.65 of a plain
# cylinder (line 4) and its base stress is 510 kPa against 780 kPa for the cylinder
# (line 5), so tapering saves stone AND lowers the load on the bottom courses. The toy
# dovetail model gives 39 kN for one block, 392 kN for a ring of ten locked together:
# ten times the resistance, simply because all ten must move at once. It is a model
# of the idea, not an engineering calculation.