Basins of attraction - another fractal

Arjen Markus (30 november 2020) I have been intrigued by the various methods to produce a fractal image for many years. While Mandelbrot and Julia sets may be the best known pictures, there is another type of such images that I had not experimented with myself: basins of attraction that occur with the Newton-Raphson method for finding the roots of an equation.

The principle is quite simple: given an equation like (z a complex number):

   z ** 3 - 1 = 0

use the Newton-Raphson method to find the roots starting from some (complex) starting point. It turns out that depending on that starting point you can end up with any of the three roots, 1, -1/2+sqrt(3)*i/2 or -1/2+sqrt(3)*i/2.

So, using this fact, you iterate a number of times and then identify which root was approached: that determines the colour for the starting point. The Newton-Rahpson method for the above equation reads:

   z(k+1) = z(k) - (z(k)**3 - 1)/ 3*z(k)**2

And the result, which tooks several minutes to compute on my computer is:

Basins of attraction - picture

The code is really straightforward:

# newton.tcl --
#     Examine the attraction basins of the Newton-Raphson iteration
#     for the equation:
#         z**3 - 1 = 0
#
package require math::complexnumbers

namespace import ::math::complexnumbers::*

# determineBasin --
#     Use the Newton-Raphson iteration and determine which root is
#     reached
#
# Arguments:
#     x, y            Cartesian coordinates of the point
#
# Result:
#     Name of the colour associated with the root
#
# Note:
#     As we do not zoom in too much, we only need a few iterations
#
#     Roots:
#     blue: 1, red: -0.5 + 0.5 sqrt(3) i, lime: -0.5 - 0.5 sqrt(3) i
#
proc determineBasin {x y} {
    set z      [complex $x  $y ]
    set a      [complex 1.0 0.0]
    set factor [complex 3.0 0.0]

    for {set i 0} {$i < 10} {incr i} {
        set z [- $z [/ [- [* $z [* $z $z]] $a] [* $factor [* $z $z]]]]
    }

    # Which basin?
    if { abs( [real $z] - 1.0 ) < 0.00001 } {
        return "#0000FF" ;# "blue"
    }

    if { abs( [real $z] + 0.5 ) < 0.00001 } {
        if { abs( [imag $z] - 0.5*sqrt(3) ) < 0.00001 } {
            return "#FF0000" ;# "red"
        } else {
            return "#00FF00" ;# "lime"
        }
    }

    return "#000000" ;# "black"
}

pack [canvas .c -width 800 -height 800]

set dxp [expr {800 / 100}]
set dyp [expr {800 / 100}]

set xmin -2.0
set xmax  2.0
set ymin -2.0
set ymax  2.0

set img [image create photo]
.c create image 0 0 -anchor nw -image $img

for {set yp 0} {$yp < 800} {incr yp} {
    for {set xp 0} {$xp < 800} {incr xp} {

        set x [expr {($xp  - 400.0) / 200.0}]
        set y [expr {(400.0 - $yp ) / 200.0}]

        if { $x == 0.0 && $y == 0.0 } {
            set x 0.00001
            set y 0.00001
        }

        set colour [determineBasin $x $y]

        $img put $colour -to $xp $yp
    }
}

The iteration procedure should be refined - there is a lot of black space left in the boundary areas. Also, zooming in on these areas should reveal the truly fractal nature: evermore complex areas where the three colours alternate ...

An interesting exercise might be to parallellise the calculation using the thread package.