function [B] = DRM(A)
 
% This program operates on a matrix A in A_n and produces
% a permutationally similar matrix that exhibits "decreasing-
% row-multiplicity", i.e., rows with higher multiplicity
% are above rows with lower multiplicity and identical rows
% are adjacent.
 
[m n] = size(A);
if (m ~= n)
    B = 0;
    'Error: matrix A must be square'
end
k = floor(sqrt(n));
if (k^2 ~= n)
    B = 0;
    'Error: matrix A must be nxn with n = k^2 for positive integer k'
end
 
uniqueRowMatrixA = zeros(n,n); % stores a copy for each unique row in A
multVectorA = zeros(1,n);      % stores number of rows equal to
corresponding
row in uniqueRowMatrixA
multMatrixA = zeros(n,k);      % stores row numbers of rows that are
equal to
the corresponding
                               %   row in uniqueRowMatrixA
counterA = 2;                  % stores current row in uniqueRowMatrixA
to
where a unique row should
                               %   be copied next
equalA = 1;                    % temp bool variable to show if two rows
are
equal up to current point
matchedRowA = 0;               % temp variable to hold row # of row in
uniqueRowMatrixA that
                               %   matches with current row of A
                               
% first row in A is first unique row
for i = 1:n
    uniqueRowMatrixA(1,i) = A(1,i);
end
multVectorA(1) = 1;
multMatrixA(1,1) = 1;
 
for i = 2:n % for each remaining row in A
 
    for r = 1:(counterA-1) % compare to each row in uniqueRowMatrixA
        for s = 1:n        % compare to element in each column
            if (A(i,s) ~= uniqueRowMatrixA(r,s))
                equalA = 0;
                break
            end
        end
        if (equalA == 1)
            matchedRowA = r;
            break
        end
        equalA = 1;
    end
    if (matchedRowA ~= 0)
        multVectorA(r) = multVectorA(r) + 1;
        multMatrixA(r,multVectorA(r)) = i;
    else
        for j = 1:n
            uniqueRowMatrixA(counterA,j) = A(i,j);
        end
        multVectorA(counterA) = multVectorA(counterA) + 1;
        multMatrixA(counterA,1) = i;
        counterA = counterA + 1;
    end
    equalA = 1;
    matchedRowA = 0;
end
 
newA = zeros(n,n);
count = 1;
count2 = 1;
 
for i = 0:(k-1)
    for j = 1:n
        if (multVectorA(j) == k - i)
            for f = 1:k
                if (multMatrixA(j,f) ~= 0)
                    listA(count) = multMatrixA(j,f);
                    count = count+1;
                end
            end
            newMultVectorA(count2) = k - i;
            newMultVectorA(count2+1) = 0;
            count2 = count2 + 1;
        end 
    end
end
 
I = eye(n,n);
P_A = I(listA,:);
 
newA = P_A * A * P_A';
B = newA;