Mandelbrot with rational exponents, without branch cut in iteration
3z[k+1] = z[k]^p + c, for p = m/n ∈ ℚ
When p ∈ ℤ, z^p is a pretty simple matter. However, with p ∉ ℤ, a problem arises: because arg(z) is ambiguous (adding any multiple of 2π yields the same number in Cartesian form), raising z to a non-integral exponent can yield multiple results (for a rational exponent, the number of distinct results is equal to the denominator n).
To deal with this, I iterate for each of n distinct arguments, between 0 and 2nπ (I only calculate the half of the image with non-negative imaginary parts, so the argument is always non-negative; the lower half is a mirror image), adding multiples of 2π. During the iteration, I consider the argument to be unambiguous and calculate purely in polar coordinates (while keeping the argument small to avoid loss of precision; adding or subtracting a multiple of 2nπ doesn't change the result in theory). This way, I can track it on the Riemann surface. And yes, it is possible to _add_ vectors in polar coordinates.
The results from the n calculations are combined, with different hues assigned to different branches, scaled for outside points according to number of iterations to escape, and all colors added. Also, I chose a denominator of 360. Because it is a highly composite number, many values get highly reduced, resulting in fewer branches that must be calculated.
Problems with the current algorithm:
* How to calculate an appropriate minimal escape radius? For integer p, 2 is fine, but for other values, it can increase quite a bit. I don't know what the relation to p is. For this sequence, I just left it at 2, so the results aren't really what I would consider completely correct.
* Currently, I calculate the colors by equally distributing different hues in HSV space, then adding colors in sRGB. This isn't really a good way to do it. Plus, the area where all branches are in the set turns out pink, instead of the fully saturated red it's supposed to be, when p ∉ ℤ.
* I'm wondering if there's a more efficient way to normalize arg(z). Currently, I just subtract 2nπ until it's less than 2nπ. Maybe dividing by 2nπ, taking the floor, and subtracting that times 2nπ?