# leidingnet.py — rekenpomp uit i11, water bij ongeveer 20 °C
import math
RHO, NU, G = 1000.0, 1.0e-6, 9.81

def factor(re):
    if re <= 0:
        return 0.0
    if re <= 2300:
        return 64 / re
    turbulent = 0.3164 * re**(-0.25)  # gladde leiding, Blasius
    if re >= 4000:
        return turbulent
    w = (re - 2300) / 1700          # alleen schatting in overgangsgebied
    return (1-w)*64/re + w*turbulent

def verlies(q, D=0.016, L=2.5, K=2.0):
    v = (q / 60000) / (math.pi*D*D/4)
    re = v*D/NU
    return (factor(re)*L/D + K)*v*v/(2*G)

def pomp(q):
    return 3.0 * (1-q/10.0)

def werkpunt(D=0.016, L=2.5, K=2.0, z=1.5):
    lo, hi = 0.0, 10.0
    for _ in range(70):
        q = (lo+hi)/2
        if pomp(q) > z + verlies(q,D,L,K):
            lo = q
        else:
            hi = q
    q = (lo+hi)/2
    return q, pomp(q)

for d in (0.008, 0.010, 0.013, 0.016):
    q,h = werkpunt(D=d)
    v = q/60000/(math.pi*d*d/4)
    print(f"D={d*1000:2.0f} mm: Q={q:.3f} L/min, H={h:.3f} m, Re={v*d/NU:.0f}")
q,h = werkpunt()
ph = RHO*G*q/60000*h
print(f"hydraulisch {ph:.3f} W, elektrisch 12 W: rendement {100*ph/12:.2f} %")
print(f"met extra filter K=20: {werkpunt(K=20)[0]:.3f} L/min")