%On teste la convergence de la méthode de newton pour "z^3-1=0" dans C
%Ce script s'inspire du script du livre de Lapresté sur Matlab, page 107
%Voir aussi http://fr.wikipedia.org/wiki/Fractale_de_Newton

clear all
[x y] = meshgrid(-.8:.0044:.6,-.79:.0044:.79); z = x+i.*y;
zs = roots([1,0,0,-1]);

%initialisation
tol=0.001; %tolérance
nmax=25;   %nombre maximum d'itérations

im = zeros(size(z)); %l'image initiale
for n=1:nmax
   test = (abs(z-zs(2))>tol) & (abs(z-zs(1))>tol) & (abs(z-zs(3))>tol); %cf. "find"
   %test(i,j) == 1 ssi la convergence n'est pas atteinte en z(i,j)
   if sum(sum(test)) < 2000
      break;
   end;
   z(test) = 2/3*z(test)  + 1./(3*z(test).^2);
   im(test) = im(test) + 1;
end

figure(1)
fact = max(im(:));
colormap gray, image((1-im/fact)*100), axis equal
axis('off');