function [s, z] = seplam_initial(A, B, prtlevel)
% search along lines parallel to the real axis passing through all 2n eigenvalues
% of A and B, looking for ALL real positive s and z on these real lines such that
% sigma_min(A - zI) = sigma_min(B - zI) = s. Return the smallest such s and
% the corresponding z. This ensures that we don't get stuck with a component of
% A's eps-pseudospectrum strictly contained inside a component of B's eps-pseudospectrum,
% for some eps, thus defeating "realimagcheck".
% could make this more efficient
% can save factor of 2 by looking along lines joining eigenvalues
% can avoid repeat checks along same line, e.g. if eigenvalues are real
eA = eig(A);
eB = eig(B);
eigtol = 1e-8;
for j=1:length(A)
[s,z] = linecheck(A, B, eA(j), eigtol);
[val,indx] = min(s);
sminA(j) = val;
zminA(j) = z(indx);
if prtlevel > 1
fprintf('\nfound initial upper bound %e on line %d for A', val, j)
end
end
for j=1:length(B)
[s,z] = linecheck(A, B, eB(j), eigtol);
[val,indx] = min(s);
sminB(j) = val;
zminB(j) = z(indx);
if prtlevel > 1
fprintf('\nfound initial upper bound %e on line %d for B', val, j)
end
end
[valA,indxA] = min(sminA);
[valB,indxB] = min(sminB);
if valA < valB
s = valA;
z = zminA(indxA);
else
s = valB;
z = zminB(indxB);
end
% possible that s is inf and z is nan, if no s and z were found on
% any of these lines