function Near2Instab(A)

[n,n] = size(A);
i = sqrt(-1);
evalues = eig(A);

% tau1(A) Discussion.
% Suppose Q'*A*Q = T is the schur decomp and the eigenvalue that is
% closest to the imaginary axis is T(k,k). Assume A is stable.
% If ek = I(:,k) then T + E1 where E1 = |Real(T(k,k)|ek*ek' has an 
% eigenvalue on the imaginary axis. It is easy to show that norm(E1,2) =
% |Real(T(k,k))| Thus
%       T + E1 = Q'*A*Q + E1  = Q'(A + E)Q  where E = Q*E1*Q'
% It follows that A+E has an eigenvalue on the imaginary axis and
% norm(E,2) = norm(E1,2) = |Real(T(k,k))| = Tau1(A)

tau1 = min(abs(real(evalues)));

m = norm(A,1);
m = 2*max(abs(evalues));   % A better range based on P2X
N = 1000;
muVals = linspace(-m,m,N);
tau2Vals = zeros(1,N);

tol = .000001*norm(A,1)/norm(A,'fro');   % For reasoning see below...
for k=1:N
    if k==1
        [U,S,V]       = jacSVD(A+muVals(1)*i*eye(n,n),tol);
    else
        Ak = A+muVals(k)*i*eye(n,n);
        % Preprocess with old U and V..
        Bk = U'*Ak*V;
        [delU,S,delV] = jacSVD(Bk,tol);
        U = U*delU; 
        V = V*delV;
    end
    tau2Vals(k) = min(abs(diag(S)));
    % 
end
[tau2,j] = min(tau2Vals);
semilogy(muVals,tau2Vals);
title(sprintf('tau1(A) = %10.2e   tau2(A) =  %10.2e   muopt = %10.5f',tau1,tau2,muVals(j)))
xlabel('mu')
ylabel('sigmamin(A+mu*i*I)')

% Suppose when JacSvD terminates the "final" A matrix has the form
% D + E where D is diagonal and E is the matrix of off diagonal elements.
% It follows that off(E) <= tol norm(A,'fro') and that D + E has the
% same singular values as A.
% Let dmin be the smallest diagonal entry of D and sigmin the smallest
% singular value of A. The theorem tells us that
%     |dmin - sigmin| <= norm(E,'fro') <= tol*norm(A,'fro')
% To guarantee that |dmin - sigmin| <= .001norm(A,1) we just insist
% that tol*norm(A,'fro') <= .001*norm(A,1). I.e., 
% set tol = .001 *norm(A,1)/norm(A,'fro').

% The .001 was a bad choice given the particular example. A better
% plot results with the .001 replaced by .000001