Python - how to compute all nth roots of a number?

Question

Is it possible to calculate n complex roots of a given number using Python? I've shortly checked it, and it looks like Python gives me wrong/incomplete answers:

(-27.0j)**(1.0/3.0) produces (2.598076211353316-1.4999999999999998j)

but proper roots should be 3 complex numbers, because every non-zero number has n different complex number nth roots. Is it possible in Python?

Answer 1

I don't think standard Python will do this unless you write a function for it, but you can do it with Numpy:

http://docs.scipy.org/doc/numpy/reference/generated/numpy.roots.html

Answer 2

There are many multi-valued complex functions - functions that can have more than one value corresponding to any point in their domain. For example: roots, logarithms, inverse trigonometric functions...

The reason these functions can have multiple values is usually because they are the inverse of a function that has multiple values in the domain map to the same value.

When doing calculations with such functions, it would be impractical to always return all possible values. For the inverse trigonometric functions, there are infinitely many possible values.

Usually the different function values can be expressed as a function of an integer parameter k. For example, the values of log z with z = r*(cos t + i*sin t is log r + i*(t + k*2*pi) with k any integer. For the nth root, it is r**(1/n)*exp(i*(t+k*2*pi)/n with k=0..n-1 inclusive.

Because returning all possible values is impractical, mathematical functions in Python and almost all other common programming languages return what's called the 'principal value' of the function. (reference) The principal value is usually the function value with k=0. Whatever choice is made, it should be stated clearly in the documentation.

So to get all the complex roots of a complex number, you just evaluate the function for all relevant values of k:

def roots(z, n):
    nthRootOfr = abs(z)**(1.0/n)
    t = phase(z)
    return map(lambda k: nthRootOfr*exp((t+2*k*pi)*1j/n), range(n))

(You'll need to import the cmath module to make this work.) This gives:

>>> roots(-27j,3)
[(2.59808-1.5j), (1.83691e-16+3j), (-2.59808-1.5j)]