
clear

f = @(x)  x.^2 - 2;
df = @(x) 2*x;

[x,k, fout, xout] = newton(f, df, 1, 1e-8, 100);




f = @(x) x.^(1/3);
df = @(x) 1/3 * x.^( -2/3);
[x,k, fout, xout] = newton(f, df, 1, 1e-8, 100);












g = @(z) z.^3 - 1;
dg = @(z) 3 * z.^2;

[x,k, fout, xout] = newton(g, dg, -1+0.2i, 1e-8, 100);










N = 100;
for j = 1:2*N+1
    for k = 1:2*N+1
        z0 = (k-N-1)/N + (j-N-1)/N * 1i;
        [z] = newton(g, dg, z0, 1e-8, 100);
        z_conv(j,k) = z;
    end
end


root_num = zeros(size(z_conv));
root_num( find(abs(z_conv-1) < 1e-7)) = 1;
root_num( find(abs(z_conv+1/2 + sqrt(3)/2 * 1i) < 1e-7)) = 2;

root_num( find(abs(z_conv+1/2 - sqrt(3)/2 * 1i) < 1e-7)) = 3;

imagesc((-N:N)/N, (-N:N)/N,root_num)
axis xy
xlabel('x')
ylabel('y')
