RAMSES	Systems Ecology
Start         Modules:   A-Z   Function   Layer        QuickRefs:   DM   AuxLib   AuxLibE   SciLib   EasyMW   MW   ISIS   RMSLib

DEFINITION MODULE Jacobi;

  (*******************************************************************

    Module  Jacobi     (Version 1.0)

      Copyright (c) 1992-2006 by Dimitrios Gyalistras and ETH Zurich.

    Purpose   Compute the eigenvalues and eigenvectors of
              a real symmetric matrix.

    Remarks   Algorithms taken from Press, W.H. et al. (1988).
              Numerical Recipes, Cambridge, pp 335 and 748-750
              (ISBN 0521 30811 9).


    Programming

      o Design
        Dimitrios Gyalistras      03/04/1992

      o Implementation
        Dimitrios Gyalistras      03/04/1992


    ETH Zurich
    Systems Ecology
    CHN E 35.1
    Universitaetstrasse 16
    8092 Zurich
    SWITZERLAND

    URLs:
        <mailto:RAMSES@env.ethz.ch>
        <http://www.sysecol.ethz.ch>
        <http://www.sysecol.ethz.ch/SimSoftware/RAMSES>


    Last revision of definition:  15/04/1992  DG

  *******************************************************************)


  CONST
    VecSize=40;

  TYPE
    Vector = ARRAY [1..VecSize] OF REAL;
    Matrix = ARRAY [1..VecSize] OF Vector;


  PROCEDURE Jacobi( VAR mat       : Matrix;
                        dim       : INTEGER;
                    VAR eigVals   : Vector;
  					VAR eigVecs   : Matrix;
  					VAR numRot    : INTEGER );
  (*
    Calculates all eigenvectors and eigenvalues of the symmetric
    dim x dim -matrix stored in 'mat'.  'eigVals' returns in its
    first 'dim' elements the eigenvalues of 'mat'.  'eigVecs' is
    a dim x dim -matrix returning the normalized eigenvectors of
    'mat' (sum of squares of vector elements = 1).  'numRot'
    returns the number of jacobi rotations which were required.

    IMPLEMENTATION RESTRICTION: After calling the procedure
    Jacobi, all elements of 'mat' above the diagonal will have
    been overwritten with 0!  No test is done against a
    non-symmetrical matrix mat, i.e. all results returned for a
    non-symmetrical mat are generally WRONG!!!
  *)

  PROCEDURE EigSort( VAR eigVals: Vector;
                     VAR eigVecs: Matrix;
                     dim: INTEGER );
  (*
    Given 'dim' eigenvalues in vector 'eigVals' and the dim x
    dim -matrix 'eigVecs' containing the corresponding
    eigenvectors (as obtained both from procedure Jacobi), this
    routine sorts the eigenvalues in descending order. The columns
    of the 'eigVecs'-matrix are rearranged correspondingly. The
    method used is straight insertion.
  *)

END Jacobi.