Simplify matrix to have ones in diagonal

I have a function:

function [A_new, b_new] = forward_elimination(A, b)
%FORWARD_ELIMINATION - Performs forward elimination to put A into unit
% upper triangular form.
%   A - original matrix of Ax = b
%   b - original vector of Ax = b
%   A_new - unit upper triangular A formed using Gaussian Elimination
%   b_new - the vector b associated with the transformed A
A_new = A;
b_new = b;
[n,n]=size(A);
if any(diag(A)==0)
  error('cannot compute')
end
for row=1:n-1
  for i=row+1:n
      factor=A(i,row)/A(row,row);
      for j=row:n
          A(i,j)=A(i,j)-factor*A(row,j);
      end
      b(i)=b(i)-factor*b(row);
  end
  A_new=A;
  b_new=b;
end
end

This function provides the correct answer, however, I want the diagonals to be one. For example, when inputting a matrix A=[1 2 3; 4 5 6; 7 8 8] and b=[1;2;3], so: [A_new,b_new]=forward_elimination(A,b)

I want to produce a simplified matrix.

In this example, A_new=[1 2 3; 0 1 2; 0 0 1] b_new=[1;2/3;0]

However, my current code produces A_new=[1 2 3; 0 -3 -6; 0 0 -1] and b_new=[1;-2;0]

How do I simplify the function even more?

Answer 1

[David Goodmanson]

Hi amintr,

You can use

d = diag(A);    
A_new = A./d;   % newer versions of Matlab with implicit expansion
b_new = d.*b;

% without implicit expansion:
n = size(A,1);
A_new = A./repmat(d,1,n);