This code implements the piping rules according to Dutch WBI standard 2017, update 2021-05-28. See this manual page
I have checked it with an existing Piping Excel sheet but feel free to test some more and give me feedback (breinbaasnl@gmail.com)
Please don't mind the Dutch / English mixups.. some words are easier in Dutch (well.. for me that is ;-)
class Piping(BaseModel):
"""Class for piping calculations
Note that this assumes a deklaag of one cohesive layer
# TODO > allow the calculation of soilstresses instead of using 1 soillayer for the deklaag
Source: Schematiseringshandleiding piping, WBI 2017, 28 mei 2021
https://www.helpdeskwater.nl/onderwerpen/waterveiligheid/primaire/beoordelen/@205752/schematiseringshandleiding-piping/
Args:
hh_polder (float): stijghoogte in de polder [m]
h_exit (float): hoogte uittredepunt maaiveld of waterpeil in sloot (denk aan de effectieve laagdikte berekening) [m]
r_exit (float): damping factor over kwellengte (0.0 = no damping, 1.0 = full damping) [-]
hh_entry (float): stijghoogte at entry point (usualy the waterlevel in the river) [m]
h_coh (float): dikte van de cohesieve laag (bovenzijde is gelijk aan h_exit) [m]
y_sat (float): saturated volumetric weight of the cohesive soillayer (will not work if there are multiple soillayers in the cohesive layer!) [kN/m3]
y_water (float): volumetric weight of water, defaults to 9.81 [kN/m3]
"""
hh_polder: float
h_exit: float
r_exit: float
hh_entry: float
h_coh: float
y_sat: float
y_water: float = 9.81
@property
def hh_exit(self) -> float:
"""Stijghoogte at exit point"""
return self.hh_polder + self.r_exit * (self.hh_entry - self.hh_polder)
def uplift(self) -> float:
"""Calculate uplift
Returns:
float: uplift factor
"""
f = (self.h_coh * (self.y_sat - self.y_water)) / self.y_water
return f - (self.hh_exit - self.h_exit)
def heave(self, i_ch: float = 0.5) -> float:
"""Calcute heave
Args:
i_ch (float): critical heave gradient, defaults to 0.5 (safe value)
Returns:
float: heave factor
"""
return i_ch - (self.hh_exit - self.h_exit) / self.h_coh
class Sellmeijer(Piping):
h_acquifer: float # dikte van de watervoerende laag [m]
k: float # Darcy permeability sand [m/s]
d70: float # 70percentiel korrelverdeling sand [m]
mp: float = 1.0 # model factor [-]
rc: float = 0.3 # reduction factor entry resistance at cohesive layer [-]
def factor_from_l(
self, l_entry_exit: float
) -> float:
"""Get the Sellmeijer factor from the distance between the entry and exit point [m]
Args:
l_entry_exit (float): Distance between entry and exit point
Returns:
float: Sellmeijer factor
"""
n = 0.25 # white coefficient [-]
rho = 37 # rolweerstandshoek sand [deg]
d70m = 2.08e-4 # d70 ref value [m]
vw = 1.33e-6 # kinematic viscosity of water [1.33e-6 (at 10C)]
g = 9.81 # gravity [m/s2]
ysp = 16.5 # vol weight sand korrels
fr = n * (ysp / self.y_water) * tan(radians(rho))
ddl = self.h_acquifer / l_entry_exit
a = pow(ddl, 2.8) - 1
b = 0.28 / a + 0.04
fg = 0.91 * pow(ddl, b)
kappa = vw / g * self.k
c = d70m / pow(kappa * l_entry_exit, 1 / 3)
d = self.d70 / d70m
fs = c * pow(d, 0.4)
dhc = l_entry_exit * fr * fs * fg
return self.mp * dhc - (self.hh_entry - self.h_exit - self.rc * self.h_coh)
Don't mind the from leveelogic stuff.. you are seeing part of a much bigger python lib for levee assessments and I copied this right out of it so you will have to do your own import but there is no dependency to any other LeveeLogic code necessary to run this code.
import pytest
from leveelogic.calculations.piping import Piping, Sellmeijer
class TestPiping:
def test_hh_exit(self):
c = Piping(
hh_polder=-1.5,
h_exit=-2.0,
r_exit=0.006,
hh_entry=0,
h_coh=2.5,
y_sat=14.0,
)
assert round(c.hh_exit, 2) == -1.49
def test_uplift(self):
c = Piping(
hh_polder=-1.5,
h_exit=-2.0,
r_exit=0.006,
hh_entry=0,
h_coh=2.5,
y_sat=14.0,
)
assert round(c.uplift(), 2) == 0.56
def test_heave(self):
c = Piping(
hh_polder=-1.5,
h_exit=-2.0,
r_exit=0.006,
hh_entry=0,
h_coh=2.5,
y_sat=14.0,
)
assert round(c.heave(), 2) == 0.30
def test_sellmeijer(self):
s = Sellmeijer(
hh_polder=-1.5,
h_exit=-2.0,
r_exit=0.006,
hh_entry=0,
h_coh=2.5,
y_sat=14.0,
h_acquifer=5,
k=4e-4,
d70=2.43e-4,
)
assert round(s.factor_from_l(l_entry_exit=9), 2) == -0.38