Your odes are generating infinite results early on. You are getting NaN for all c values.
WIth the below small modification to break early when NaN is detected, infinite values show up for c
isfirstnan = true;
i=1;
error=0.0001;
C=zeros(100,1);
aC=zeros(100,1);
for a=0.01:0.01:1
span=3.6;
yspan=[-span span];
k=0;
u0=[0.09,a*0.09];
c0=complex(0,0.1);
sol=ode45(@(y,u) func(y,u,c0,a),yspan,u0);
U=deval(sol,span);
A0=(U(2)/U(1))+a;
c1=complex(0,0.12);
sol=ode45(@(y,u) func(y,u,c1,a),yspan,u0);
U=deval(sol,span);
A1=(U(2)/U(1))+a;
c=(((c0*A1)-(c1*A0)/(A1-A0)));
while k<1000
sol=ode45(@(y,u) func(y,u,c,a),yspan,u0);
U=deval(sol,span);
A0=A1;
A1=(U(2)/U(1))+a;
c0=c1;
c1=c;
errR=abs(real(c1-c0));
errI=abs(imag(c1-c0));
if((errR<error)&&(errI<error))
break;
end
c=(((c0*A1)-(c1*A0)/(A1-A0)));
if (any(~isfinite(real(c))) || any(~isfinite(imag(c))))
if isfirstnan
fprintf('got first non-finite for c at a = %g, k = %g, i = %g\n', a, k, i);
disp(c);
isfirstnan = false;
end
break
end
k=k+1;
end
c;
k;
C(i)=imag(c);
aC(i)=a*imag(c);
i=i+1;
end
whos C, min(C), max(C)
whos aC, min(aC), max(aC)
%{
t=linspace(0.01,1,100)';
plot(t,C)
hold on
plot(t,aC)
hold off
ylim([0 1.1]);
xlim([0.01 1]);
legend('Ci','αCi');
%}
function dudy=func(y,u,c,a)
dudy=zeros(2,1);
dudy(1,1)=u(2);
dudy(2,1)=(a*a+2*tanh(y)*(tanh(y)*tanh(y)-1)/(tanh(y)-c))*u(1);
end
