Двойная лексикография c сортировка бинарных (0-1) матриц

Question

Я хочу отсортировать двоичную матрицу так, чтобы ее столбцы и строки были в лексикографическом порядке путем переключения строк и столбцов. Другими словами, мне нужно двойное лексическое упорядочение матрицы, состоящей из нулей и единиц, в то время как записи всех строк и столбцов остаются в соответствующей строке / столбце. Это эффективно заполняет нижнюю правую часть результирующей матрицы как можно большим числом 1 с, не разбивая столбцы / строки. (В моем случае мне нужно заполнить верхнюю левую сторону, поэтому я просто перевернул результат.)

Может ли кто-нибудь дать мне ссылку на код, делающий это? Очень ценным бонусом было бы то, что алгоритм работает с тензорами любого измерения.

Сначала я попытался отсортировать все лексикографически все строки, а затем все столбцы, повторяя это до тех пор, пока полученная матрица больше не изменится. Я сравнил это с первой сортировкой столбцов и использовал лучший результат. Однако при этом результат, казалось, застрял в каком-то локальном оптимуме и не дал наилучшего возможного результата. Поэтому я искал существующие алгоритмы, но не нашел ничего реализованного и быстрого для использования.

Я нашел старую статью на sciencedirect.com [1], которая делает именно это в O (mn) но есть проблемы с его реализацией. Может быть, кто-то с доступом может помочь мне с кодом, особенно с

• pointers, которые предназначены для отслеживания блоков столбцов, которые еще не были разделены при вызове dlrefine.
• И с заполнением bucket, где он может несколько раз потерпеть неудачу.

[1]: Spinrad, Jeremy P. "Вдвойне лексическое упорядочение плотных матриц 0–1 «. Письма обработки информации 45.5 (1993): 229-235.

Вот ошибочный код Matlab из моей попытки реализации:

function [A,I]=dlpartition(M)
    Rout={};
    RP={1:size(M,1)};
    Rlast=1;
    CP={1:size(M,2)};
    pointers=ones(1,size(M,1));
    while Rlast>0
        R=RP{Rlast};
        C=[];
        if any(pointers(R)>0)
            Clast=pointers(R);
            Clast=min(Clast(Clast>0));
            C=CP{Clast}; %the last column block unexamined by R
        end
        if isempty(C)
            Rlast=Rlast-1;
            Rout=[R Rout];
        else
            [R,C,pointers]=dlrefine(M,R,C,pointers-Clast+1);
            pointers=pointers+Clast-1;
            RP=[RP{1:Rlast-1} R RP{Rlast+1:end}];
            Rlast=Rlast-1+length(R);
            CP=[CP{1:Clast-1} C CP{Clast+1:end}];
        end
    end

    %prepare output to match my needs
    I=reshape(1:numel(M),size(M));
    I=I(flip([Rout{:}]),:);
    I=I(:,flip([CP{:}]));
    A=M(I);
end

function [R,C, pointers]=dlrefine(M,R,C,pointers)
%REFINE Spinrad 1992
%M        Binary matrix to sort
%R        Array of rows in current block
%C        Array of columns in current block
%pointers Array of pointers from all rows to the
%         column block that wasn't split yet by it
%         with offset, such that 1 is the current
%         block C
    C={C};
    for r=R
        c=length(C);
        C2=[];
        if pointers(r) >= c %if Cc has not been split by r
            while isempty(C2) && c > 0
                C1=intersect(C{c},find(M(r,:)));
                C2=setdiff(C{c},C1);
                if ~isempty(C1) && ~isempty(C2)
                    C={C{1:c-1} C2 C1 C{c+1:end}};
                    pointers=pointers+(pointers>=c);
                    %pointers(r)=c-1; %TESTING (normally uncommented)
                end
                %if isempty(C2) %TESTING (normally uncommented)
                    c=c-1;
                %end %TESTING (normally uncommented)
            end
            pointers(r)=c; %TESTING (normally commented out)
        end
    end
    bucket=cell(1,length(C));
    failed=[];
    for r=R
        success=false;
        for sub=length(C):-1:1
            if any(ismember(C{sub},find(1-M(r,:))))
                bucket{sub}=[bucket{sub} r];
                success=true;
                break;
            end
        end
        if ~success
            failed=[failed r];
        end
    end
    R=bucket(end:-1:1);
    R=[R(~cellfun('isempty',R)) failed];
end

Для удобства вы можете запустить dlpartition функция с той же матрицей, что и в статье:

M=[0 1 0 1 1 0;
   1 0 1 1 1 0;
   0 1 1 1 1 0;
   1 1 0 0 0 1;
   0 1 0 1 0 1;
   1 0 0 1 1 1];
[A I]=dlpartition(M)

Редактировать:

Код работает для небольшого примера, приведенного в статье, и выдает:

A =  1 1 1 1 0 0
     1 1 1 0 0 0
     1 1 0 1 1 0
     1 1 0 0 1 1
     1 0 1 0 0 1
     0 0 1 0 1 1

Но он не может генерировать глобальный оптимум для произвольных матриц. Например:

M=[1 0 1 0 0
   1 1 0 0 1
   1 0 0 0 1
   0 0 1 1 1
   0 0 1 0 1];
[A I]=dlpartition(M)
%this gives            but should give
%A = 1 1 0 0 0         1 1 1 0 0
%    1 0 1 1 0         1 1 0 0 0
%    1 0 1 0 0         1 0 0 1 0
%    0 1 1 0 1         0 1 0 1 1
%    0 1 1 0 0         0 1 0 1 0