terminate loop diverging to inf
显示 更早的评论
I am trying to iterate a basins of attraction map, where the values diverge to either 0 or inf. When they diverge to inf the loop is continuous. I would like to terminate it at a particularly large value and store the index value as p(k,kk)=2 and then continue iterations, and similarly when it diverges to 0 to store as p(k,kk)=1
Any help is massively appreciated !
NT = 1024;
N = 1024;
r = 0.9;
omega = (sqrt(5)+1)/2;
phi = 2*pi*omega;
eps = 1e-6;
x0_vec = linspace(-1.0,1.8,N); y0_vec = linspace(-1.0,1.8,N);
z0_vec = complex(x0_vec, y0_vec);
P = zeros(N,N);
for k = 1:length(y0_vec)
y0 = y0_vec(k);
for kk = 1:length(x0_vec)
x0 = x0_vec(kk);
z = z0_vec;
for j = 1:NT
z(j+1) = z(j)^2 + r*exp(1i*phi)*z(j);
%stop loop continuing to inf
if abs(z(j+1)) <= eps
disp('Winning attractor: at 0');
P(k,kk) = 1;
else
disp('Winning attractor: at infinity');
P(k,kk) = 2;
end
end
end
end
1 个评论
dpb
2023-11-24
MAGICNUMBEER=99999999999999999999;
....
%loop iterator here...
if newvalue>=MAGICNUMBER
P(k,kk)=2
break % breaks innermost loop only...
end
采纳的回答
I'm not sure, but I guess you meant something like this:
NT = 1024;
N = 2048;
r = 0.9;
omega = (sqrt(5)+1)/2;
phi = 2*pi*omega;
eps = 1e-6;
x0_vec = linspace(-1.0,1.8,N); y0_vec = linspace(-1.0,1.8,N);
P = zeros(N,N);
for i = 1:numel(x0_vec)
for j = 1:numel(y0_vec)
z0 = complex(x0_vec(i),y0_vec(j));
for k = 1:NT
z0=z0^2+r*exp(1i*phi)*z0;
if abs(z0) <= eps
%disp('Winning attractor: at 0');
P(i,j)=1;
break
end
if abs(z0) > 1e8
%disp('Winning attractor: at infinity');
P(i,j)=2;
break
end
end
end
end
contourf(x0_vec,y0_vec,P)
colorbar
sum(P(:)==1)
sum(P(:)==2)
sum(P(:)==1)+sum(P(:)==2)
numel(x0_vec)*numel(y0_vec)
更多回答（1 个）
There is one critical point is the eps (MATLAB built in var ==> not to use) = 1e-6 is not met. The smallest abs value of z is 0.012; therefore, EPS_lim (different name instead of eps) = 0.0012; Why not ot use vectorization instead of the loop. The [for end] loop is very slow and just displays the phrases in disp() for over million times, e.g.:
NT = 1024;
N = 1024;
r = 0.9;
omega = (sqrt(5)+1)/2;
phi = 2*pi*omega;
EPS_lim = 0.012;
x0_vec = linspace(-1.0,1.8,N);
y0_vec = linspace(-1.0,1.8,N);
z0_vec = complex(x0_vec, y0_vec);
P = zeros(N,N);
z = z0_vec;
z(2:end+1) = z(1:end).^2 + r*exp(1i*phi)*z(1:end);
IDX1 = abs(z(2:end)) <= EPS_lim;
P(:,IDX1)=1;
IDX2 = P==0;
P(IDX2)=2;
disp('Number of Winning attractor: at 0');
sum((P(:)==1))
disp('Winning attractor: at infinity');
sum((P(:)==2))
1 个评论
Walter Roberson
2023-11-24
By the way, matlab eps is a function rather than a variable. You can pass a parameter to it, and the result is the epsilon relative to the input value.
类别
在 帮助中心 和 File Exchange 中查找有关 Image Processing Toolbox 的更多信息
Translated by 
