function gerschgorin(A)
%GERSCHGORIN(A)
%   Plots the Gerschgorin circles of a matrix A,
%   both based on row-sums (red) as on column-sums (blue),
%   together with the eigenvalues of the matrix.
%
%   (C) 2010/12/01 by Benjamin Seibold

lambda = eig(A); % compute eigenvalues
t = linspace(0,2*pi,100); cx = cos(t); cy = sin(t); % circle points
n = size(A,1); name = 'rows'; clf
for k = 1:2 % rows/columns
    subplot(2,1,k), hold on
    for i = 1:n % loop over all rows/columns
        c = A(i,i); % center of circle
        r(i) = sum(abs(A(i,[1:i-1,i+1:end]))); % radius of circle
        h(1) = patch(c+r(i)*cx,r(i)*cy,[k-1 0 2-k],'edgecolor','none');
    end
    re = max(r)*.05; % radius of eigenvalue markers
    for i = 1:n % loop over all rows/columns
        h(2) = patch(real(lambda(i))+re*cx,imag(lambda(i))+re*cy,...
        [0 0 0],'edgecolor','none');
    end
    hold off, axis equal, ax{k} = axis;
    legend(h,['circles from ',name],'eigenvalues')
    name = 'columns'; A = A';
end
for k = 1:2 % make common axes
    subplot(2,1,k)
    axis([min(ax{1}(1),ax{2}(1)),max(ax{1}(2),ax{2}(2)),...
          min(ax{1}(3),ax{2}(3)),max(ax{1}(4),ax{2}(4))])
end

