Newton Basins MATLAB Code
Here we give only the code for construction of the Newton basins of a cubicpolynomial.
For a general polynomial of given degree its value and derivativederiv can be entered explicitly in the for loopthat carries out the iteration of Newton's method.
function w = newtbas(poly,wind,colors,count)
% Uses count iterations to draw Newton basins
% in the viewing window wind for the polynomial
% with descending vector poly of coefficients,
% using the RGB color matrix colors. Try
%
% poly = [1 -3 0 1];
% wind = [-2 4 -2.25 2.25];
% colors =
% [0.0 0.0 0.0; % black
% 0.9 0.9 0.3; % yellow
% 0.8 0.0 0.1; % red
% 0.2 0.7 0.3]; % green
% count = 50;
%
% for the cubic polynomial x^3 - 3x^2 + 1.
r = roots(poly); % Find polynomial roots
s = size(r); n = s(1);
r1 = r(1); r2 = r(2); r3 = r(3);
a = wind(1); b = wind(2);
c = wind(3); d = wind(4);
k = 300; % no of gridpoints
x = a : (b-a)/(k-1) : b;
y = c : (d-c)/(k-1) : d;
[x,y] = meshgrid(x,y);
z = x + i*y;
clear x y
for n = 1:count
s1 = z - r1; s2 = z - r2; s3 = z - r3;
% value = s1.*s2.*s3;
% deriv = s2.*s3 + s1.*s3 + s1.*s2;
% z = z - value./deriv; % Newton's iteration
z = z - s1.*s2.*s3./(s2.*s3 + s1.*s3 + s1.*s2);
end
clear value deriv s1 s2 s3
w1 = (abs(z - r1) < 0.01); % Determine which root
w2 = (abs(z - r2) < 0.01); % the initial guess z
w3 = (abs(z - r3) < 0.01); % yields upon iteration
w = 1 + w1 + 2*w2 + 3*w3; % convert roots to colors
image(w), colormap(newtmap), axis off