def laminar_entry_thermal_Hausen(Re=None, Pr=None, L=None, Di=None):
    r'''Calculates average internal convection Nusselt number for laminar flows
    in pipe during the thermal entry region according to [1]_ as shown in
    [2]_ and cited by [3]_.

    .. math::
        Nu_D=3.66+\frac{0.0668\frac{D}{L}Re_{D}Pr}{1+0.04{(\frac{D}{L}
        Re_{D}Pr)}^{2/3}}

    Parameters
    ----------
    Re : float
        Reynolds number, [-]
    Pr : float
        Prandtl number, [-]
    L : float
        Length of pipe [m]
    Di : float
        Diameter of pipe [m]

    Returns
    -------
    Nu : float
        Nusselt number, [-]

    Notes
    -----
    If Pr >> 1, (5 is a common requirement) this equation also applies to flows
    with developing velocity profile.
    As L gets larger, this equation  becomes the constant-temperature Nusselt
    number.

    Examples
    --------
    >>> laminar_entry_thermal_Hausen(Re=100000, Pr=1.1, L=5, Di=.5)
    39.01352358988535

    References
    ----------
    .. [1] Hausen, H. Darstellung des Warmeuberganges in Rohren durch
       verallgeminerte Potenzbeziehungen, Z. Ver deutsch. Ing Beih.
       Verfahrenstech., 4, 91-98, 1943
    .. [2] W. M. Kays. 1953. Numerical Solutions for Laminar Flow Heat Transfer
       in Circular Tubes.
    .. [3] Bergman, Theodore L., Adrienne S. Lavine, Frank P. Incropera, and
       David P. DeWitt. Introduction to Heat Transfer. 6E.
       Hoboken, NJ: Wiley, 2011.
    '''
    Gz = Di/L*Re*Pr
    return 3.66 + (0.0668*Gz)/(1+0.04*(Gz)**(2/3.))