Mandelbrot Set and Newton's Method Basins of Attraction

September 19 2005 Alec Mihailovs 4100
Maple
2
This blog entry evolved from my reply to Moira Chas's post. I want to thank her for initiating such an interesting topic. Usually Mandelbrot set is drawn in Maple using plot3d command. That also can be done using densityplot. In the example below I use mandelbrot from John Oprea's worksheet,
mandelbrot := proc(x, y) 
local c, z, m; 
c := evalf(x+y*I); z := c; 
for m to 50 while abs(z) < 2 do z := z^2+c od; 
m end;

plots[densityplot](mandelbrot,-2 .. 0.7, -1.35 .. 1.35,
s_tyle=patchnogrid,colorstyle=HUE,numpoints=62500,axes=none);
One needs to replace s_tyle with style to make this working (there is a problem with using 'style' inside <pre>.) Changing 50 in the mandelbrot procedure to 30, 40, or other number, produces slightly different pictures. The picture also looks interesting in RGB,
plots[densityplot](mandelbrot,-2 .. 0.7, -1.35 .. 1.35,
s_tyle=patchnogrid,colorstyle=RGB,numpoints=62500,axes=none);
Newton's Method Basins of Attraction can be implemented in a similar way, using the following procedure,
NF := proc(f, n, r1, r2, num, cs)
local F, DF, L, d, Newt, x, k, L1;
    F := unapply(f, op(indets(f, name)));
    DF := D(F);
    L := [fsolve(f, op(indets(f, name)), complex)];
    d := nops(L);
    Newt := proc(x0)
            x := evalf(x0);
            to n do x := evalf(x - F(x)/DF(x)) end do;
            L1 := map(y -> abs(y - x), L);
            member(min(op(L1)), L1, 'k');
            evalf((k - 1)/d)
        end proc;
    print(plots[densityplot]('Newt'(a + b*I), a = r1, b = r2,
        s_tyle = patchnogrid, colorstyle = cs, numpoints = num,
        axes = none))
end proc:
The procedure is not very efficient and it takes some time to produce pictures. I hope to improve it some time (especially if somebody would make some good suggestions). Here are few pictures obtained using this procedure,
NF(z^2-1,5,-2..2,-2..2,625,HUE);
NF(z^3-1,5,-2..2,-2..2,62500,HUE);
NF(z^4-1,6,-2..2,-2..2,62500,HUE);
NF(z^5-1,6,-2..2,-2..2,62500,HUE);
NF(z^6-1,6,-2..2,-2..2,62500,HUE);
NF(x^7-2*x^6+3*x^4-4,6,-2..2,-2..2,62500,HUE);

Please Wait...