MATLAB Answers

Reduced row echelon form technique

97 views (last 30 days)
I am trying to use a code to calculate the reduced row echelon form of a matrix without the function rref. I make a random matrix A and and then make a matrix new_A = (A-lambda*I). The intent is to eventually find the nullspace of new_A without the null function. When comparing the code to the rref command the results match a majority of the time, however on certain occasions there is a discrepancy in the final column. My understanding is that rref must be unique for a given matrix. Any ideas on what could be the source?
clear all
clc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% MAKES RANDOM MATRIX A
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
m = 3;
casenum = randi(3,1)
if casenum == 1
A = randn(m);
end
if casenum == 2
A = i*randn(m);
end
if casenum == 3
pownum = randi(2,m);
A = randn(m) + randn(m).*(i.^pownum);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% MAKES MATRIX new_A FROM (A-lambda*I) & CALCULATES MATLAB NULL
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
I = eye(m,m);
mat_eigs = eig(A);
new_A = A - mat_eigs(1)*I % currently only uses 1st eigenvalue to test
mat_null = null(new_A);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FINDS RREF OF new_A
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
mat_rref = rref(new_A)
tol = 1e-12;
i = 1;
j = 1;
while (i <= m) && (j <= m)
% Find value and index of largest element in the remainder of column j
[p,k] = max(abs(new_A(i:m,j)));
k = k+i-1;
if (p <= tol)
% The column is negligible, zero it out
new_A(i:m,j) = 0;
j = j + 1;
else
% Swap i-th and k-th rows
new_A([i k],j:m) = new_A([k i],j:m);
% Divide the pivot row by the pivot element
Ai = new_A(i,j:m) / new_A(i,j);
% Subtract multiples of the pivot row from all the other rows
new_A(:,j:m) = new_A(:,j:m) - new_A(:,j)*Ai;
new_A(i,j:m) = Ai;
i = i + 1;
j = j + 1;
end
end
code_rref = new_A

  0 Comments

Sign in to comment.

Accepted Answer

Richard Brown
Richard Brown on 7 Aug 2013
It'll be a difference between yours and MATLAB's choice of tolerance. Your matrix should be rank 2, however I notice that MATLAB is often computing the RREF to be the identity.
If you read the docs for rref they specify what they use as tolerance:
max(size(A))*eps *norm(A,inf))

  1 Comment

Lingling Fan
Lingling Fan on 13 Jun 2019
Hi,
Thanks for your answer!
How could we change the tolerance of rref, as I did not find in docs.
Is there something like "options" like fmincon setup for tolerance?
Thanks in advance!

Sign in to comment.

More Answers (0)

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!