User:DX-MON/Complex Number Newton-Raphson differentiation

This article outlines how to do Complex Number based Newton-Raphson iteration, and differentiating the complex number function to be iterated.

The Newton-Raphson iteration method is defined as: $f(z)=z-{\frac {g(z)}{g'(z)}}$
for any arbitrary function $g(z)\,$

I am defining the following for this differentiation:

• $g(z)=z^{3}-1\,$
  
• $z=x+yi\,$

Expanding the function g(z)

The logical expansion of g(z) is:
$g(z)=g(x,yi)=x^{3}+3x^{2}yi-3xy^{2}-y^{3}i-1\,$

which can be partially differentiated with respect to x or y as below:

${\frac {\partial {g}}{\partial {x}}}=3x^{2}+6xyi-3y^{2}$
${\frac {\partial {g}}{\partial {y}}}=3x^{2}i-6xy-3y^{2}i$

${\frac {\partial {g}}{\partial {x}}}$ can be used as the real component of the output of the function $g'(z)\,$ and ${\frac {\partial {g}}{\partial {y}}}$ the imaginary component of the output, except this is only true after switching two terms of the equations so that one is purely in terms of y and x and the other in terms of yi and x. This yealds the following for ${\frac {\partial {g}}{\partial {x}}}$ and ${\frac {\partial {g}}{\partial {y}}}$:

${\frac {\partial {g}}{\partial {x}}}=3x^{2}-6xy-3y^{2}$
${\frac {\partial {g}}{\partial {y}}}=3x^{2}i+6xyi-3y^{2}i$

Defining g'(z)

With the two partial differentiations, we can say that $g'(z)\,$ should be defined as:
$g'(z)=g'(x,y)=(3x^{2}-6xy-3y^{2},3x^{2}i+6xyi-3y^{2}i)\,$
where the output is represented by a vector who's x coordinate is the real value of the complex number and the y coordinate is the imaginary one in terms of $i\,$. If only life were quite so simple as this though. The combination of the two partial differentiations is only two thirds of the picture - to have a working equation for programming this, we need to perform the full differentiation too (differentiate with respect to x and y together), so: ${\frac {dg}{dxy}}=3x^{2}+6xi-6y-3y^{2}i$

putting this into our definition of $g'(x,y)\,$, we get: $g'(z)=g'(x,y)=(3x^{2}-6xy-3y^{2}+3x^{2}-6y,3x^{2}i+6xyi-3y^{2}i+6yi-3y^{2}i)\,$

Conclusion

now we have a full definition of both $g(z)\,$ and $g'(z)\,$, we can perform the iteration, code for doing so is presented below in Python:

from sys import argv
from PIL import Image, ImageDraw
from math import sqrt
from cmath import sin

def g_norm(z):
    return ((z ** 3) - 1)

def g_diff(z):
    return complex((3 * (z.real ** 2)) - (6 * (z.real * z.imag)) - (3 * (z.imag ** 2)),
        (3 * (z.real ** 2)) + (6 * (z.real * z.imag)) - (2 * (z.imag ** 2)))

def Julia(width, height, scale):
    Set_img = Image.new("RGB", (width, height), (255, 255, 255))
    Set = ImageDraw.Draw(Set_img)
    maxiters = 50

    w = width / 2
    h = height / 2

    print "Calculating Newton-Raphson Julia Set and drawing/colouring image"
    for x in range(-w, w):
        for y in range(-h, h):
            z = complex(x, y)
            z *= scale
            i = 0
            while i < maxiters:
                i += 1
                lz = z
                try:
                    z = z - (g_norm(z) / g_diff(z))
                except:
                    break
                if abs(z) == abs(lz):
                    break
            i = maxiters - i
            colour = int((float(i) / float(maxiters)) * 255.0)
            if colour > 10:
                if z.imag == 0.0:
                    colour = (0, colour, 0)
                elif z.imag > 0.0:
                    colour = (0, 0, colour)
                else:
                    colour = (colour, 0, 0)
            else:
                colour = (colour, colour, colour)
            Set.point((x + w, y + h), fill=colour)
            if y == (h - 1):
                print (float(x + w) / float(width)) * 100.0, "% done"
    print "100% done"
    del Set
    Set_img.save("Julia.png", "PNG")

if  __name__ == '__main__':
    if len(argv) > 3:
        Julia(int(argv[1]), int(argv[2]), float(argv[3]))
    else:
        Julia(200, 200, 1.0)

There is one problem with this though, if run it will not produce the expected Julia Set and Fatou Set. An approximation of the correct set is obtained when g_diff(z) is defined as

def g_diff(z):
    return (3 * (z ** 2))