Efficient Resolution of the Colebrook Equation

colebrook.py

from math import log

def vel(q, a):
    """
    Mean Velocity (q[m3/s], a[m2])
    q is flow rate
    a is cross-sectional area
    """
    return q / a

def Re(vel, Dh, nu=0.000001):
    """
    Reynolds number (vel[m/s], Dh[m], nu=0,000001[m2/s])
    vel is mean velocity
    Dh is hydraulic diameter
    """
    return vel * Dh / nu

def k(roughness, Dh):
    """
    Relative Roughness (roughness[mm], Dh[m])    
    """
    return roughness / (Dh * 1000)


def Colebrook(Re, k):
    """
    Re is the Reynolds number
    k is relative roughness (k = e/Dh)
    http://math.unice.fr/~didierc/DidPublis/ICR_2009.pdf
    """
    X1 = k * Re * 0.123968186335418
    X2 = log(Re) - 0.779397488455682
    F = X2 - 0.2
    E = (log(X1 + F) + F - X2) / (1 + X1 + F)
    F = F - (1 + X1 + F + 0.5 * E) * E * (X1 + F) / (1 + X1 + F + E * (1 + E / 3))
    E = (log(X1 + F) + F - X2) / (1 + X1 + F)
    F = F - (1 + X1 + F + 0.5 * E) * E * (X1 + F) / (1 + X1 + F + E * (1 + E / 3))
    F = 1.15129254649702 / F
    return F * F