% FIGURE 3.3.4
      function EIG = DEGNON(A,N)
%
%  FUNCTION DEGNON SOLVES THE EIGENVALUE PROBLEM
%
%                A*X = LAMBDA*X
%
%    WHERE A IS A GENERAL REAL MATRIX.
%
%
%  ARGUMENTS
%
%             ON INPUT                          ON OUTPUT
%             --------                          ---------
%
%    A      - THE N BY N MATRIX.
%
%    N      - THE SIZE OF MATRIX A.
%
%    EIG    -                                   A COMPLEX N-VECTOR
%                                               CONTAINING THE EIGEN-
%                                               VALUES OF A.
%
%-----------------------------------------------------------------------
%                              EPS = MACHINE FLOATING POINT RELATIVE
%                                    PRECISION
% *****************************
      EPS = eps;
% *****************************
%                              AMAX = MAXIMUM ELEMENT OF A
      AMAX = 0.0;
      for I=1:N
      for J=1:N
         AMAX = max(AMAX,abs(A(I,J)));
      end
      end
      ERRLIM = sqrt(EPS)*AMAX;
%                              REDUCTION TO HESSENBERG FORM
      A = HESSQ(A,N);
%                              REDUCTION TO QUASI-TRIANGULAR FORM
      A = QR(A,N,ERRLIM);
%                              EXTRACT EIGENVALUES OF QUASI-TRIANGULAR
%                              MATRIX
      I = 1;
      while (I <= N-1)
         if (A(I+1,I) == 0.0)
%                              1 BY 1 BLOCK ON DIAGONAL
            EIG(I) = A(I,I);
            I = I+1;
         else
%                              2 BY 2 BLOCK ON DIAGONAL
            DISC = (A(I,I)-A(I+1,I+1))^2 + 4.0*A(I,I+1)*A(I+1,I);
            TERM =  0.5*(A(I,I)+A(I+1,I+1));
            if (DISC >= 0.0)
               EIG(I)  = TERM + 0.5*sqrt(DISC);
               EIG(I+1)= TERM - 0.5*sqrt(DISC);
            else
               EIG(I)  = TERM + 0.5*sqrt(-DISC)*i;
               EIG(I+1)= TERM - 0.5*sqrt(-DISC)*i;
            end
            I = I+2;
         end
      end
      if (I == N)
         EIG(N) = A(N,N);
      end