Hi all, @c.poupin
I intend to design the foundation of an elevated water tank, consisting of a raft with the shape shown in the image below. To do so, I need to determine the raft diameter Dr by solving the following inequality:
𝑁𝑇∙𝐷𝑟 ≥ 8∙𝑀 (where NT and M are inputs in my code)
Based on my hand calculations, after taking into account all the necessary parameters to evaluate NT, the final equation to be solved is a cubic polynomial with the following expression:
5694·Dr³ - 0.942·Dr² +2540.505·Dr - 73612.480 = 0
The real solution obtained from my hand calculation is: Dr = 17,36 m
However, when I run the code below, the result is:
cubic polynomial equation: 6.676·Dr³ - 0.942·Dr² +2540.505·Dr - 73612.480 = 0
reel solution: Dr = 16.744 m
What could be causing this discrepancy, and how can it be fixed?
Raft shaoe:
Please check my code here:
import sys
import os
sys.path.append(r'C:\Users\nono\AppData\Local\Programs\Python\Python38\Lib\site-packages')
import numpy as np
import math
def vertical_load(H, a, b, c, di, de, gb, gs, P0, Dr):
# calculating the entire vertical force at the base
# H is embedment depth
# de: shaft outer diameter
# di: shaft inner diameter
# gb is the concrete density = 2.5T/m3
# gs is the soil density = 1.9T/m3
# P0 vertical load from the superstructure
# Dr raft diameter to find
h = H / 4
d1 = de + 2*a
d2 = d1 + 2*b
d3 = Dr - 2*c
# raft area calculation
s1 = de**2 - di**2
s2 = Dr**2
s3 = (d3**2 - d2**2) / 2
s4 = d2**2
s5 = (d2**2 - d1**2) / 2
s6 = d1**2
# raft weight calculation
Pr = gb * h * math.pi / 4 * (s1 + s2 + s3 + s4 + s5 + s6)
# soil area calculation
t1 = Dr**2 - d3**2
t2 = (d3**2 - d2**2) / 2
t3 = Dr**2 - d2**2
t4 = (d2**2 - d1**2) / 2
t5 = Dr**2 - de**2
# soil weight
Pt = gs * h * math.pi / 4 * (t1 + t2 + t3 + t4 + t5)
# Total vertical load
NT = P0 + Pr + Pt
return NT
def Equi_inequality(Dr, H, a, b, c, di, de, gb, gs, P0, M):
NT = vertical_load(H, a, b, c, di, de, gb, gs, P0, Dr)
return NT * Dr - 8 * M
def polynom(H, a, b, c, di, de, gb, gs, P0, M):
Dr_vals = np.array([0.0, 1.0, 2.0, 3.0])
F_vals = np.array([
Equi_inequality(Dr, H, a, b, c, di, de, gb, gs, P0, M)
for Dr in Dr_vals
])
A = np.array([
[0, 0, 0, 1],
[1, 1, 1, 1],
[8, 4, 2, 1],
[27, 9, 3, 1]
])
coeffs = np.linalg.solve(A, F_vals)
return coeffs
def solve_Dr(H, a, b, c, di, de, gb, gs, P0, M):
a3, a2, a1, a0 = polynom(
H, a, b, c, di, de, gb, gs, P0, M
)
roots = np.roots([a3, a2, a1, a0])
real_roots = [
r.real for r in roots
if abs(r.imag) < 1e-8 and r.real > 0
]
if not real_roots:
raise ValueError("no solution found")
return min(real_roots), (a3, a2, a1, a0)
# Inputs
H = 4.00
a = 1.05
b = 2.00
c = 1.00
di = 6.80
de = 7.60
gb = 2.5
gs = 1.9
P0 = 2492.52
M = 9201.56
Dr, coeffs = solve_Dr(H, a, b, c, di, de, gb, gs, P0, M)
print(f"{coeffs[0]:.3f}·Dr³ {coeffs[1]:+.3f}·Dr² "
f"{coeffs[2]:+.3f}·Dr {coeffs[3]:+.3f}")
print(f"Dr = {Dr:.3f} m")
OUT = Dr, coeffs = solve_Dr(H, a, b, c, di, de, gb, gs, P0, M)
Thanks.
