home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l460 / 3.ddi / DEMOS.DI$ / RREFMOVI.M < prev    next >
Encoding:
Text File  |  1993-03-07  |  2.2 KB  |  75 lines

  1. function ratmovie(A,tol)
  2. %RREFMOVIE Movie of the computation of the reduced row echelon form.
  3. %    RREFMOVIE(A) produces the reduced row echelon form of A.
  4. %    RREFMOVIE, by itself, supplies its own 8-by-6 matrix with rank 4.
  5. %    RREFMOVIE(A,tol) uses the given tolerance in the rank tests.
  6.  
  7. %    Copyright (c) 1984-93 by The MathWorks, Inc.
  8.  
  9. % Sample matrix if none was provided.
  10. if nargin < 1
  11.  
  12.     A = [ 9     4     1     6    12     7
  13.           4     0     4    15     1    14
  14.           7     0     7     8    10     9
  15.          16     0    16     3    13     2
  16.           0     2    -4     0     0     0
  17.           0     6   -12     0     0     0
  18.           9     0     9     6    12     7
  19.           5     0     5    10     8    11];
  20. end
  21.  
  22. format rat
  23. more off
  24. clc
  25. home
  26. disp('  Original matrix')
  27. A
  28. disp('Press any key to continue. . .'), pause
  29. [m,n] = size(A);
  30.  
  31. % Compute the default tolerance if none was provided.
  32. if (nargin < 2), tol = max([m,n])*eps*norm(A,'inf'); end
  33.  
  34. % Loop over the entire matrix.
  35. i = 1;
  36. j = 1;
  37. k = 0;
  38. while (i <= m) & (j <= n)
  39.    % Find value and index of largest element in the remainder of column j.
  40.    [p,k] = max(abs(A(i:m,j))); k = k+i-1;
  41.    if (p <= tol)
  42.       % The column is negligible, zero it out.
  43.       home
  44.       disp(['  column ' int2str(j) ' is negligible'])
  45.       A(i:m,j) = zeros(m-i+1,1)
  46.       disp('Press any key to continue. . .'), pause
  47.       j = j + 1;
  48.    else
  49.       if i ~= k
  50.          % Swap i-th and k-th rows.
  51.          home
  52.          disp(['  swap rows ' int2str(i) ' and ' int2str(k) blanks(10)])
  53.          A([i k],:) = A([k i],:)
  54.          disp('Press any key to continue. . .'), pause
  55.       end
  56.       % Divide the pivot row by the pivot element.
  57.       home
  58.       disp(['  pivot = A(' int2str(i) ',' int2str(j) ')' blanks(10)])
  59.       A(i,j:n) = A(i,j:n)/A(i,j)
  60.       % Subtract multiples of the pivot row from all the other rows.
  61.       home
  62.       disp(['  eliminate in column ' int2str(j) blanks(10)])
  63.       for k = 1:m
  64.          if k ~= i
  65.             home
  66.             disp(' ')
  67.             A(k,j:n) = A(k,j:n) - A(k,j)*A(i,j:n)
  68.          end
  69.       end
  70.       disp('Press any key to continue. . .'), pause
  71.       i = i + 1;
  72.       j = j + 1;
  73.    end
  74. end
  75.