Skip to content

Instantly share code, notes, and snippets.

@breinbaas
Last active September 1, 2022 16:31
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save breinbaas/af9a10c150ec72639ca58fe786186244 to your computer and use it in GitHub Desktop.
Save breinbaas/af9a10c150ec72639ca58fe786186244 to your computer and use it in GitHub Desktop.

Piping code

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)

The code (scroll down to see some test samples)

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)

Tests of use cases

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
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment