clear all;
%choose a set of values of k between 0 and 1 to test.
ks = [0.1:0.1:0.9]';

allpoles = []; %we'll store all the poles we will compute in this array
for k_index=1:length(ks)
    
    k = ks(k_index);
    ki_min = (4 + 3*k - 2*sqrt(2*k - k^2))/k;
    ki_max = (4 + 3*k + 2*sqrt(2*k - k^2))/k;
    kis = linspace(ki_min+1e-3,ki_max -1e-3,9);
    
    for ki_index = 1:length(kis)
        ki = kis(ki_index);
        poles = roots([1,4+2*k, k*ki*2]);
        allpoles = cat(1,allpoles,poles(:));
    end;
end;

%check if the |p - c|<1 condition is satisfied:

c = -3+ 2*j;
index_of_poles_positive_imaginary_part = find(imag(allpoles)>0);
p_s = allpoles(index_of_poles_positive_imaginary_part);

distances = abs(p_s - c);
disp(['maximum value of |p-c| is ',num2str(max(distances))]);