def RD(a,b,c):
    """
    Compute Carlson's elliptic integral RD(a, b, c).

                   3  /\oo              dx
    RD(a, b, c) =  -  |   ------------------------------
                   2  |          1/2       1/2       3/2
                     \/0  (x + a)   (x + b)   (x + c)

    The parameters a, b, and c may be complex numbers.
    """
    A, B, C = a, b, c
    sum = 0.0
    for k in range(5):
        gms = (A**0.5)*(B**0.5) + (A**0.5)*(C**0.5) + (B**0.5)*(C**0.5)
        sum += 1 / (4**k * C**0.5 * (C+gms))
        A, B, C = (A+gms)/4, (B+gms)/4, (C+gms)/4
    avg0, avg = (a+b+3*c)/5., (A+B+3*C)/5
    X = (avg0 - a)/(avg * 4**5)
    Y = (avg0 - b)/(avg * 4**5)
    Z = (-X-Y)/3
    E2, E3, E4, E5 = X*Y-6*Z*Z, (3*X*Y-8*Z*Z)*Z, 3*(X*Y-Z*Z)*Z*Z, X*Y*Z*Z*Z
    return 3*sum + avg**-1.5/4**5 * (1 + E2*(-3/14. + E2*9/88. - E3*9/52.)
                                       + E3/6. - E4*3/22. + E5*3/26.)