function [V,H] = qriterS(A);
%
%
%   Performs QR iteration in real arithmetic
%
%   Input:  A -- a real n by n symmetric matrix
%
%   Output: V -- a real n by n orthogonal matrix
%
%           H -- a real n by n  diagonal matrix
%   
%                AV = VH
%
%           Note:  In real life, you would take advantage of symmetry.
%                  This code does not do that.
% 
%   D.C. Sorensen
%   2 Mar 2000
%

%
    [n,n] = size(A);
%
%   Reduce A to tridiagonal form with full Arnoldi process
%
    v = randn(n,1);
    [V,H] = ArnoldiC(A,n,v);

    k1 = 1; k2 = n;
    while(k2 > 1)
       if (k2 == 2),
%
%         Diagonalize the leading 2 by 2 block after all
%         the rest of the eigenvalues have converged.
%
          [Q,D] = eig(H(1:2,1:2));

          H(1:2,:) = Q'*H(1:2,:);
          H(:,1:2) = H(:,1:2)*Q;
          V(:,1:2) = V(:,1:2)*Q;
          H(2,1) = 0;
          k2 = k2 - 1;
       else
%
%      Select the Wilkinson shift 
%
          mu = eig(H(k2-1:k2,k2-1:k2));
          if (abs(mu(1) - H(k2,k2)) < abs(mu(2) - H(k2,k2)) ),
             [V,H] = qrstep(V,H,mu(1),k1,k2);
          else
             [V,H] = qrstep(V,H,mu(2),k1,k2);
          end
%
%      Deflate: Set off-diagonal to zero if it is small compared to
%               the two adjacent diagonal elements

          if (abs(H(k2,k2-1)) < 10*eps*max(abs(H(k2-1,k2-1)),abs(H(k2,k2)))),
             H(k2,k2-1) = 0;
             H(k2-1,k2) = 0;
             k2 = k2 - 1;
          end
       end

%
%      ---------------Graphics Display of Convergence--------------------
%


       x = ones(n*(n+1)/2 , 1); y = x;
       phd = zeros(n-1,1);

       m = 16;
       colors = bluemap(m);
       axis([0 n+1 0 n+1]);
       x = [-1  -1  n+2 n+2];
       y = [-1  n+2 -1  n+2];
       ph1 = plot(x,y,'o');
       axis tight;
       set(ph1,'MarkeredgeColor',[1 1 1])
       hold on;
       indx = 1;
       Hmax = max(abs(diag(H,-1)));
       Hmax = max(Hmax,1);
       x = [1:n];
       y = [n:-1:1];
       for j = 1:n-1, 
           jj = min(n-1,j);
           i = min(n,jj+1);
           gval = max(eps, abs(H(i,jj)/Hmax));
           gval = -log(gval);
           gval = max(1,ceil(gval));
           gval = min(16,gval);
           gval = 16 - gval + 1;
           ii = n-i+1; jj = j; 

           phd(j) = plot([jj,jj+1],[ii, ii+1],'o');
           axis tight;
           set(phd(j),'Color',colors(m,:))
           set(phd(j),'MarkerSize',9)
           set(phd(j),'MarkeredgeColor',[1 1 1])

           if(abs(H(i,j)) > 0),
              set(phd(j),'MarkerFaceColor', colors(gval,:));
           end
       end
       ph = plot(x,y,'o');
       axis tight;
       set(ph,'Color',colors(m,:))
       set(ph,'MarkerSize',9)
       set(ph,'MarkeredgeColor',[1 1 1])
       set(ph,'MarkerFaceColor',colors(m,:))
       pause(.1)
%      pause
       hold off
       axis tight;
    end