Ответ 1

(Обновленные результаты тестирования, примеры прогона и фрагменты кода ниже.)

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

Рассмотрим пример матрицы 7x7; в начале состояние:

0,0,0,0,0,0,0  

означает, что имеется семь смежных неиспользуемых столбцов. После добавления двух в первую строку состояние может быть, например:

0,1,0,0,1,0,0  

с двумя столбцами, которые теперь имеют один в них. После добавления ко второй строке состояние может быть, например:

0,1,1,0,1,0,1  

После заполнения трех строк существует вероятность того, что столбец будет иметь максимум два; это эффективно разбивает матрицу на две независимые зоны:

1,1,1,0,2,0,1  ->  1,1,1,0 + 0,1  

Эти зоны независимы в том смысле, что правило без смежных не действует при добавлении в разные зоны, а порядок зон не влияет на количество решений.

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

2,1,1,0 + 0,1  

где 2 означает, что этот столбец использовался в текущей строке (а не в двух столбцах). На следующем шаге мы должны преобразовать двойки обратно в единицы.

Кроме того, мы можем также отразить отдельные группы, чтобы поместить их в свои лексикографически наименьшие обозначения:

2,1,1,0 + 0,1  ->  0,1,1,2 + 0,1  

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

0,0 + 0,1 + 0,0,2 + 0,1,0 + 0,1,0,1  

Затем при вычислении количества решений, возникающих в каждом состоянии, мы можем использовать memoization, используя каноническое обозначение каждого состояния в качестве ключа.

Создание словаря состояний и количества решений для каждого из них нужно сделать только один раз, а таблица для более крупных матриц, вероятно, может быть использована и для меньших матриц.

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

Обратите внимание, что каждое состояние и соответствующий ключ могут встречаться только в определенной строке, поэтому вы можете хранить ключи в отдельных словарях в строке.

РЕЗУЛЬТАТЫ ИСПЫТАНИЙ

Первый тест с использованием неоптимизированного JavaScript дал очень многообещающие результаты. При динамическом программировании вычисление количества решений для матрицы 10x10 теперь занимает второе место, где алгоритм грубой силы занял несколько часов (и это часть алгоритма, который нужно выполнить только один раз). Размер словаря с сигнатурами и числами решений растет с уменьшающимся фактором, приближающимся к 2,5 для каждого шага по размеру; время его генерации растет с коэффициентом около 3.

Это количество решений, состояний, подписей (общий размер словарей) и максимальное количество подписей на строку (наибольший словарь на строку), которые создаются:

size                           unique solutions          states    signatures    max/row

 4x4                                               2            9           6           2  
 5x5                                              16           73          26           8  
 6x6                                             722          514         107          40  
 7x7                                          33,988        2,870         411         152  
 8x8                                       2,215,764       13,485       1,411         596  
 9x9                                     179,431,924       56,375       4,510       1,983  
10x10                                 17,849,077,140      218,038      13,453       5,672  
11x11                              2,138,979,146,276      801,266      38,314      14,491  
12x12                            304,243,884,374,412    2,847,885     104,764      35,803  
13x13                         50,702,643,217,809,908    9,901,431     278,561      96,414  
14x14                      9,789,567,606,147,948,364   33,911,578     723,306     238,359  
15x15                  2,168,538,331,223,656,364,084  114,897,838   1,845,861     548,409  
16x16                546,386,962,452,256,865,969,596                4,952,501   1,444,487  
17x17            155,420,047,516,794,379,573,558,433               12,837,870   3,754,040  
18x18         48,614,566,676,379,251,956,711,945,475               31,452,747   8,992,972  
19x19     17,139,174,923,928,277,182,879,888,254,495               74,818,773  20,929,008  
20x20  6,688,262,914,418,168,812,086,412,204,858,650              175,678,000  50,094,203

(Дополнительные результаты, полученные с помощью С++, с использованием простой 128-битной целочисленной реализации.)

ПРИМЕР

Словарь для матрицы 5x5 выглядит следующим образом:

row 0:  00000  -> 16        row 3:  101    ->  0
                                    1112   ->  1
row 1:  20002  ->  2                1121   ->  1
        00202  ->  4                1+01   ->  0
        02002  ->  2                11+12  ->  2
        02020  ->  2                1+121  ->  1
                                    0+1+1  ->  0
row 2:  10212  ->  1                1+112  ->  1
        12012  ->  1
        12021  ->  2        row 4:  0      ->  0
        12102  ->  1                11     ->  0
        21012  ->  0                12     ->  0
        02121  ->  3                1+1    ->  1
        01212  ->  1                1+2    ->  0

Общее число решений составляет 16; если мы произвольно выбираем число от 0 до 15, например. 13, мы можем найти соответствующее (т.е. 14-е) решение следующим образом:

state:      00000  
options:    10100  10010  10001  01010  01001  00101  
signature:  00202  02002  20002  02020  02002  00202  
solutions:    4      2      2      2      2      4  

Это говорит нам о том, что 14-е решение является вторым решением опции 00101. Следующий шаг:

state:      00101  
options:    10010  01010  
signature:  12102  02121  
solutions:    1      3  

Это говорит нам о том, что второе решение является первым решением опции 01010. Следующий шаг:

state:      01111  
options:    10100  10001  00101  
signature:  11+12  1112   1+01  
solutions:    2      1      0  

Это говорит нам о том, что 1-е решение является 1-м решением опции 10100. Следующий шаг:

state:      11211  
options:    01010  01001  
signature:  1+1    1+1  
solutions:    1      1  

Это говорит нам о том, что 1-е решения являются 1-м решением опции 01010. Последний шаг:

state:      12221  
options:    10001  

И матрица 5x5, соответствующая произвольно выбранному числу 13, равна:

0 0 1 0 1  
0 1 0 1 0  
1 0 1 0 0
0 1 0 1 0  
1 0 0 0 1  

И вот пример quick'n'dirty code; запустить фрагмент для создания словаря подписи и словаря решений и создать случайную матрицу 10x10 (для создания словаря требуется секунда, а после этого он генерирует случайные решения за полмиллиона секунд):

function signature(state, prev) {
    var zones = [], zone = [];
    for (var i = 0; i < state.length; i++) {
        if (state[i] == 2) {
            if (zone.length) zones.push(mirror(zone));
            zone = [];
        }
        else if (prev[i]) zone.push(3);
        else zone.push(state[i]);
    }
    if (zone.length) zones.push(mirror(zone));
    zones.sort(function(a,b) {return a.length - b.length || a - b;});
    return zones.length ? zones.join("2") : "2";

    function mirror(zone) {
        var ltr = zone.join('');
        zone.reverse();
        var rtl = zone.join('');
        return (ltr < rtl) ? ltr : rtl;
    }
}

function memoize(n) {
    var memo = [], empty = [];
    for (var i = 0; i <= n; i++) memo[i] = [];
    for (var i = 0; i < n; i++) empty[i] = 0;
    memo[0][signature(empty, empty)] = next_row(empty, empty, 1);
    return memo;

    function next_row(state, prev, row) {
        if (row > n) return 1;
        var solutions = 0;
        for (var i = 0; i < n - 2; i++) {
            if (state[i] == 2 || prev[i] == 1) continue;
            for (var j = i + 2; j < n; j++) {
                if (state[j] == 2 || prev[j] == 1) continue;
                var s = state.slice(), p = empty.slice();
                ++s[i]; ++s[j]; ++p[i]; ++p[j];
                var sig = signature(s, p);
                var sol = memo[row][sig];
                if (sol == undefined) 
                    memo[row][sig] = sol = next_row(s, p, row + 1);
                solutions += sol;
            }
        }
        return solutions;
    }
}

function random_matrix(n, memo) {
    var matrix = [], empty = [], state = [], prev = [];
    for (var i = 0; i < n; i++) empty[i] = state[i] = prev[i] = 0;
    var total = memo[0][signature(empty, empty)];
    var pick = Math.floor(Math.random() * total);
    document.write("solution " + pick.toLocaleString('en-US') + 
        " from a total of " + total.toLocaleString('en-US') + "<br>");
    for (var row = 1; row <= n; row++) {
        var options = find_options(state, prev);
        for (var i in options) {
            var state_copy = state.slice();
            for (var j in state_copy) state_copy[j] += options[i][j];
            var sig = signature(state_copy, options[i]);
            var solutions = memo[row][sig];
            if (pick < solutions) {
                matrix.push(options[i].slice());
                prev = options[i].slice();
                state = state_copy.slice();
                break;
            }
            else pick -= solutions;
        }
    }
    return matrix;

    function find_options(state, prev) {
        var options = [];
        for (var i = 0; i < n - 2; i++) {
            if (state[i] == 2 || prev[i] == 1) continue;
            for (var j = i + 2; j < n; j++) {
                if (state[j] == 2 || prev[j] == 1) continue;
                var option = empty.slice();
                ++option[i]; ++option[j];
                options.push(option);
            }
        }
        return options;
    }
}

var size = 10;
var memo = memoize(size);
var matrix = random_matrix(size, memo);
for (var row in matrix) document.write(matrix[row] + "<br>");

Ответ 2

Введение

Вот какой прототип-подход, пытаясь решить более общую задачу единая комбинаторная выборка, которая для нашего подхода здесь означает: мы можем использовать этот подход для всего, что мы можем сформулировать как SAT-problem.

Он не использует вашу проблему напрямую и делает тяжелый объезд. Этот объезд для SAT-проблемы может помочь в отношении теории (более мощные общие теоретические результаты) и эффективности (SAT-решатели).

Говоря это, это не подход, если вы хотите пробовать в течение нескольких секунд или меньше (в моих экспериментах), по крайней мере, будучи обеспокоены единообразием.

Теория

Подход, основанный на результатах теории сложности, следует эта работа:

GOMES, Carla P.; САБАРВАЛ, Ашиш; SELMAN, Барт. Почти однородная выборка комбинаторных пространств с использованием ограничений XOR. В: Достижения в нейронных системах обработки информации. 2007. С. 481-488.

Основная идея:

• сформулировать проблему как проблему SAT
• добавить произвольно сгенерированные xors к задаче (действуя только на решающие переменные, что важно на практике)
  • это уменьшит количество решений (некоторые решения станут невозможными)
  • сделайте это в цикле (с настроенными параметрами), пока не останется только одно решение!
    • поиск какого-либо решения выполняется с помощью SAT-решателей или # SAT-solvers (= подсчет модели)
    • Если есть более одного решения: никакие xors не будут добавлены, но будет выполнен полный перезапуск: добавьте random-xors к стартовой проблеме!

Гарантии:

• при правильной настройке параметров этот подход позволяет получить почти однородную выборку
  • эта настройка может быть дорогостоящей, поскольку она основана на аппроксимации числа возможных решений
  • эмпирически это также может быть дорогостоящим!
  • Ответ Ante, в котором упоминается числовая последовательность A001499, фактически дает хорошую верхнюю границу пространства решений (поскольку она просто игнорирует ограничения смежности!)

Недостатки:

• неэффективен для больших проблем (в общем, не обязательно по сравнению с альтернативами, такими как MCMC и co.)
  • необходимо изменить/уменьшить параметры для создания образцов
  • те уменьшенные параметры теряют теоретические гарантии
  • но эмпирически: хорошие результаты все еще возможны!

Параметры:

На практике параметры:

• N: количество добавленных xors
• L: минимальное количество переменных - часть одного xor-ограничения
• U: максимальное количество переменных, принадлежащих одному xor-ограничению

N важно уменьшить количество возможных решений. Учитывая, что N константа, другие переменные, конечно, также влияют на это.

Теория говорит (если я правильно интерпретирую), что мы должны использовать L = R = 0.5 * # dec-vars.

Это невозможно на практике здесь, поскольку xor-ограничения сильно вредят SAT-решателям!

Здесь еще несколько научных слайдов о влиянии L и U.

Они называют xors размером 8-20 short-XORS, в то время как нам нужно будет использовать еще более короткие версии позже

Реализация

Окончательная версия

Вот довольно хакерская реализация в python, используя скрипты XorSample из здесь.

В основе используемого SAT-решателя лежит Cryptominisat.

Код в основном сводится к:

• Превратите проблему в конъюнктивную нормальную форму
  • как DIMACS-CNF
• Реализовать выборку:
  • Вызывает XorSample (основанный на трубах + файл)
  • Вызов SAT-решателя (на основе файлов)
• Добавить образцы в некоторый файл для последующего анализа

Код: (надеюсь, я уже предупреждал вас о качестве кода)

from itertools import count
from time import time
import subprocess
import numpy as np
import os
import shelve
import uuid
import pickle
from random import SystemRandom
cryptogen = SystemRandom()

""" Helper functions """
# K-ARY CONSTRAINT GENERATION
# ###########################
# SINZ, Carsten. Towards an optimal CNF encoding of boolean cardinality constraints.
# CP, 2005, 3709. Jg., S. 827-831.

def next_var_index(start):
    next_var = start
    while(True):
        yield next_var
        next_var += 1

class s_index():
    def __init__(self, start_index):
        self.firstEnvVar = start_index

    def next(self,i,j,k):
        return self.firstEnvVar + i*k +j

def gen_seq_circuit(k, input_indices, next_var_index_gen):
    cnf_string = ''
    s_index_gen = s_index(next_var_index_gen.next())

    # write clauses of first partial sum (i.e. i=0)
    cnf_string += (str(-input_indices[0]) + ' ' + str(s_index_gen.next(0,0,k)) + ' 0\n')
    for i in range(1, k):
        cnf_string += (str(-s_index_gen.next(0, i, k)) + ' 0\n')

    # write clauses for general case (i.e. 0 < i < n-1)
    for i in range(1, len(input_indices)-1):
        cnf_string += (str(-input_indices[i]) + ' ' + str(s_index_gen.next(i, 0, k)) + ' 0\n')
        cnf_string += (str(-s_index_gen.next(i-1, 0, k)) + ' ' + str(s_index_gen.next(i, 0, k)) + ' 0\n')
        for u in range(1, k):
            cnf_string += (str(-input_indices[i]) + ' ' + str(-s_index_gen.next(i-1, u-1, k)) + ' ' + str(s_index_gen.next(i, u, k)) + ' 0\n')
            cnf_string += (str(-s_index_gen.next(i-1, u, k)) + ' ' + str(s_index_gen.next(i, u, k)) + ' 0\n')
        cnf_string += (str(-input_indices[i]) + ' ' + str(-s_index_gen.next(i-1, k-1, k)) + ' 0\n')

    # last clause for last variable
    cnf_string += (str(-input_indices[-1]) + ' ' + str(-s_index_gen.next(len(input_indices)-2, k-1, k)) + ' 0\n')

    return (cnf_string, (len(input_indices)-1)*k, 2*len(input_indices)*k + len(input_indices) - 3*k - 1)

# K=2 clause GENERATION
# #####################
def gen_at_most_2_constraints(vars, start_var):
    constraint_string = ''
    used_clauses = 0
    used_vars = 0
    index_gen = next_var_index(start_var)
    circuit = gen_seq_circuit(2, vars, index_gen)
    constraint_string += circuit[0]
    used_clauses += circuit[2]
    used_vars += circuit[1]
    start_var += circuit[1]

    return [constraint_string, used_clauses, used_vars, start_var]

def gen_at_least_2_constraints(vars, start_var):
    k = len(vars) - 2
    vars = [-var for var in vars]

    constraint_string = ''
    used_clauses = 0
    used_vars = 0
    index_gen = next_var_index(start_var)
    circuit = gen_seq_circuit(k, vars, index_gen)
    constraint_string += circuit[0]
    used_clauses += circuit[2]
    used_vars += circuit[1]
    start_var += circuit[1]

    return [constraint_string, used_clauses, used_vars, start_var]

# Adjacency conflicts
# ###################
def get_all_adjacency_conflicts_4_neighborhood(N, X):
    conflicts = set()
    for x in range(N):
        for y in range(N):
            if x < (N-1):
                conflicts.add(((x,y),(x+1,y)))
            if y < (N-1):
                conflicts.add(((x,y),(x,y+1)))

    cnf = ''  # slow string appends
    for (var_a, var_b) in conflicts:
        var_a_ = X[var_a]
        var_b_ = X[var_b]
        cnf += '-' + var_a_ + ' ' + '-' + var_b_ + ' 0 \n'

    return cnf, len(conflicts)

# Build SAT-CNF
  #############
def build_cnf(N, verbose=False):
    var_counter = count(1)
    N_CLAUSES = 0
    X = np.zeros((N, N), dtype=object)
    for a in range(N):
        for b in range(N):
            X[a,b] = str(next(var_counter))

    # Adjacency constraints
    CNF, N_CLAUSES = get_all_adjacency_conflicts_4_neighborhood(N, X)

    # k=2 constraints
    NEXT_VAR = N*N+1

    for row in range(N):
        constraint_string, used_clauses, used_vars, NEXT_VAR = gen_at_most_2_constraints(X[row, :].astype(int).tolist(), NEXT_VAR)
        N_CLAUSES += used_clauses
        CNF += constraint_string

        constraint_string, used_clauses, used_vars, NEXT_VAR = gen_at_least_2_constraints(X[row, :].astype(int).tolist(), NEXT_VAR)
        N_CLAUSES += used_clauses
        CNF += constraint_string

    for col in range(N):
        constraint_string, used_clauses, used_vars, NEXT_VAR = gen_at_most_2_constraints(X[:, col].astype(int).tolist(), NEXT_VAR)
        N_CLAUSES += used_clauses
        CNF += constraint_string

        constraint_string, used_clauses, used_vars, NEXT_VAR = gen_at_least_2_constraints(X[:, col].astype(int).tolist(), NEXT_VAR)
        N_CLAUSES += used_clauses
        CNF += constraint_string

    # build final cnf
    CNF = 'p cnf ' + str(NEXT_VAR-1) + ' ' + str(N_CLAUSES) + '\n' + CNF

    return X, CNF, NEXT_VAR-1


# External tools
# ##############
def get_random_xor_problem(CNF_IN_fp, N_DEC_VARS, N_ALL_VARS, s, min_l, max_l):
    # .cnf not part of arg!
    p = subprocess.Popen(['./gen-wff', CNF_IN_fp,
                          str(N_DEC_VARS), str(N_ALL_VARS),
                          str(s), str(min_l), str(max_l), 'xored'],
                          stdin=subprocess.PIPE, stdout=subprocess.PIPE, stderr=subprocess.PIPE)
    result = p.communicate()

    os.remove(CNF_IN_fp + '-str-xored.xor')  # file not needed
    return CNF_IN_fp + '-str-xored.cnf'

def solve(CNF_IN_fp, N_DEC_VARS):
    seed = cryptogen.randint(0, 2147483647)  # actually no reason to do it; but can't hurt either
    p = subprocess.Popen(["./cryptominisat5", '-t', '4', '-r', str(seed), CNF_IN_fp], stdin=subprocess.PIPE, stdout=subprocess.PIPE)
    result = p.communicate()[0]

    sat_line = result.find( SATISFIABLE')

    if sat_line != -1:
        # solution found!
        vars = parse_solution(result)[:N_DEC_VARS]

        # forbid solution (DeMorgan)
        negated_vars = list(map(lambda x: x*(-1), vars))
        with open(CNF_IN_fp, 'a') as f:
            f.write( (str(negated_vars)[1:-1] + ' 0\n').replace(',', ''))

        # assume solve is treating last constraint despite not changing header!
        # solve again

        seed = cryptogen.randint(0, 2147483647)
        p = subprocess.Popen(["./cryptominisat5", '-t', '4', '-r', str(seed), CNF_IN_fp], stdin=subprocess.PIPE, stdout=subprocess.PIPE)
        result = p.communicate()[0]
        sat_line = result.find( SATISFIABLE')
        if sat_line != -1:
            os.remove(CNF_IN_fp)  # not needed anymore
            return True, False, None
        else:
            return True, True, vars
    else:
        return False, False, None

def parse_solution(output):
    # assumes there is one
    vars = []
    for line in output.split("\n"):
        if line:
            if line[0] == 'v':
                line_vars = list(map(lambda x: int(x), line.split()[1:]))
                vars.extend(line_vars)
    return vars

# Core-algorithm
# ##############
def xorsample(X, CNF_IN_fp, N_DEC_VARS, N_VARS, s, min_l, max_l):
    start_time = time()
    while True:
        # add s random XOR constraints to F
        xored_cnf_fp = get_random_xor_problem(CNF_IN_fp, N_DEC_VARS, N_VARS, s, min_l, max_l)
        state_lvl1, state_lvl2, var_sol = solve(xored_cnf_fp, N_DEC_VARS)

        print('------------')

        if state_lvl1 and state_lvl2:
            print('FOUND')

            d = shelve.open('N_15_70_4_6_TO_PLOT')
            d[str(uuid.uuid4())] = (pickle.dumps(var_sol), time() - start_time)
            d.close()

            return True

        else:
            if state_lvl1:
                print('sol not unique')
            else:
                print('no sol found')

        print('------------')

""" Run """
N = 15
N_DEC_VARS = N*N
X, CNF, N_VARS = build_cnf(N)

with open('my_problem.cnf', 'w') as f:
    f.write(CNF)

counter = 0
while True:
    print('sample: ', counter)
    xorsample(X, 'my_problem', N_DEC_VARS, N_VARS, 70, 4, 6)
    counter += 1

Результат будет выглядеть (удалены некоторые предупреждения):

------------
no sol found
------------
------------
no sol found
------------
------------
no sol found
------------
------------
sol not unique
------------
------------
FOUND

Ядро: CNF-формулировка

Введем одну переменную для каждой ячейки матрицы. N = 20 означает 400 двоичных переменных.

Adjancency:

Предварительно расчитайте все конфликты, уменьшающие симметрию, и добавьте конфликтные предложения.

Основная теория:

a -> !b
<->
!a v !b (propositional logic)

Строка/Кольцевая мощность:

Это сложно выразить в CNF, а наивные подходы нуждаются в экспоненциальном числе ограничений.

Мы используем некоторую кодировку на основе сумматора (SINZ, Carsten. К оптимальной кодировке CNF ограничений булевой мощности), которая вводит новые вспомогательные -variables.

Примечание:

sum(var_set) <= k
<->
sum(negated(var_set)) >= len(var_set) - k

Эти SAT-кодировки могут быть помещены в точные счетчики моделей (для малых N, например, 9). Число решений равно Ante, что является убедительным доказательством правильного преобразования!

Есть также интересные приблизительные модели-счетчики (также сильно основанные на xor-ограничениях), такие как approxMC, который показывает еще одну вещь, которую мы может делать с SAT-формулировкой. Но на практике я не смог их использовать (approxMC = autoconf; no comment).

Другие эксперименты

Я также создал версию, используя pblib, чтобы использовать более мощные формулировки мощности для формулировки SAT-CNF. Я не пытался использовать API на С++, но только уменьшенный pbencoder, который автоматически выбирает лучшую кодировку, которая была намного хуже, чем моя кодировка, используемая выше (что лучше всего по-прежнему является проблемой исследования, часто даже избыточными ограничениями может помочь).

Эмпирический анализ

Чтобы получить некоторый размер выборки (учитывая мое терпение), я только рассчитал образцы для N = 15. В этом случае мы использовали:

• N = 70 xors
• L, U = 4,6

Я также вычислил некоторые образцы для N = 20 с (100,3,6), но это занимает несколько минут, и мы уменьшили нижнюю границу!

Визуализация

Вот некоторые анимации (усиление отношений любви-ненависти с matplotlib):

Изменить: И (уменьшенное) сравнение с единообразной выборкой грубой силы с N = 5 (NXOR, L, U = 4, 10, 30):

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

Теория

Статистический анализ, вероятно, трудно сделать, поскольку основная проблема имеет такой комбинаторный характер. Это даже не совсем очевидно, как должен выглядеть последний файл-PDF. В случае N = нечетный, он, вероятно, неравномерен и выглядит как шахматная доска (я наблюдал за грубой силой N = 5).

Мы можем быть уверены в (imho): симметрии!

Учитывая матрицу cell-PDF, мы должны ожидать, что матрица симметрична (A = A.T). Это проверяется при визуализации, и накладывается эвклидово-норма различий во времени.

Мы можем сделать то же самое на другом наблюдении: наблюдаемые пары.

При N = 3 мы можем наблюдать следующие пары:

• 0,1
• 0,2
• 1,2

Теперь мы можем делать это для каждой строки и столбца и также должны ожидать симметрии!

К сожалению, вероятно, нелегко сказать что-то о дисперсии и, следовательно, о необходимых образцах, чтобы говорить о доверии!

Наблюдение

Согласно моему упрощенному восприятию, текущие образцы и cell-PDF выглядят хорошо, хотя конвергенция еще не достигнута (или мы далеки от единообразия).

Более важным аспектом, вероятно, являются две нормы, красиво уменьшающиеся до 0. (да, можно было бы настроить какой-то алгоритм для этого путем переноса с помощью prob = 0.5, но это не делается здесь, так как это победит его цель).

Потенциальные следующие шаги

• Параметры настройки
• Проверьте подход, используя # SAT-solvers/Model-counters вместо SAT-solvers
• Попробуйте различные формулировки CNF, особенно в отношении кодирования мощности и xor-encodings
  • XorSample по умолчанию использует tseitin-like encoding, чтобы обойти экспоненциально растущий
    • для меньших xors (как используется) может быть хорошей идеей использовать наивное кодирование (которое распространяется быстрее)
      • XorSample поддерживает это теоретически; но script работают по-разному на практике
      • Cryptominisat известен благодаря выделенной XOR-обработке (поскольку он был создан для анализа криптографии, включая много xors) и мог бы получить что-то наивным кодированием (поскольку вывод xors из взорванных CNF намного сложнее).
• Более статистический анализ
• Избавьтесь от скриптов XorSample (shell + perl...)

Резюме

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

Ответ 3

Всего несколько мыслей. Число матриц, удовлетворяющих условиям для n <= 10:

3  0
4  2
5  16
6  722
7  33988
8  2215764
9  179431924
10 17849077140

К сожалению, нет последовательности с этими числами в OEIS.

Существует одно подобное (A001499), без условий для соседних. Количество nxn-матриц в этом случае является "порядка" как число A001499 (n-1) x (n-1) матриц. Этого можно ожидать, так как число способов заполнения одной строки в этом случае, позиция 2 в n местах с по крайней мере одним нолем между ними ((n-1) выберите 2). То же, что и в позиции 2 в (n-1) местах без ограничений.

Я не думаю, что между этими матрицами порядка n и матрицей A001499 порядка n-1 существует простая связь, что означает, что если у нас есть матрица A001499, мы можем построить некоторые из этих матриц.

При этом при n = 20 число матриц составляет > 10 ^ 30. Довольно много: -/

Ответ 4

Это решение использует рекурсию для того, чтобы установить ячейку матрицы по одной. Если случайное блуждание заканчивается невозможным решением, то мы откатываемся на один шаг в дереве и продолжаем случайное блуждание.

Алгоритм эффективен, и я думаю, что сгенерированные данные очень равновероятны.

package rndsqmatrix;

import java.util.ArrayList;
import java.util.Collections;
import java.util.List;
import java.util.stream.IntStream;

public class RndSqMatrix {

    /**
     * Generate a random matrix
     * @param size the size of the matrix
     * @return the matrix encoded in 1d array i=(x+y*size)
     */
    public static int[] generate(final int size) {
        return generate(size, new int[size * size], new int[size],
                new int[size]);
    }

    /**
     * Build a matrix recursivly with a random walk 
     * @param size the size of the matrix
     * @param matrix the matrix encoded in 1d array i=(x+y*size)
     * @param rowSum
     * @param colSum
     * @return 
     */
    private static int[] generate(final int size, final int[] matrix,
            final int[] rowSum, final int[] colSum) {

        // generate list of valid positions
        final List<Integer> positions = new ArrayList();
        for (int y = 0; y < size; y++) {
            if (rowSum[y] < 2) {
                for (int x = 0; x < size; x++) {
                    if (colSum[x] < 2) {
                        final int p = x + y * size;
                        if (matrix[p] == 0
                                && (x == 0 || matrix[p - 1] == 0)
                                && (x == size - 1 || matrix[p + 1] == 0)
                                && (y == 0 || matrix[p - size] == 0)
                                && (y == size - 1 || matrix[p + size] == 0)) {
                            positions.add(p);
                        }
                    }
                }
            }
        }

        // no valid positions ?
        if (positions.isEmpty()) {

            // if the matrix is incomplete => return null
            for (int i = 0; i < size; i++) {
                if (rowSum[i] != 2 || colSum[i] != 2) {
                    return null;
                }
            }

            // the matrix is complete => return it
            return matrix;
        }

        // random walk 
        Collections.shuffle(positions);
        for (int p : positions) {
            // set '1' and continue recursivly the exploration
            matrix[p] = 1;
            rowSum[p / size]++;
            colSum[p % size]++;
            final int[] solMatrix = generate(size, matrix, rowSum, colSum);
            if (solMatrix != null) {
                return solMatrix;
            }

            // rollback 
            matrix[p] = 0;
            rowSum[p / size]--;
            colSum[p % size]--;
        }

        // we can't find a valid matrix from here => return null
        return null;
    }

    public static void printMatrix(final int size, final int[] matrix) {
        for (int y = 0; y < size; y++) {
            for (int x = 0; x < size; x++) {
                System.out.print(matrix[x + y * size]);
                System.out.print(" ");
            }
            System.out.println();
        }
    }

    public static void printStatistics(final int size, final int count) {
        final int sumMatrix[] = new int[size * size];
        for (int i = 0; i < count; i++) {
            final int[] matrix = generate(size);
            for (int j = 0; j < sumMatrix.length; j++) {
                sumMatrix[j] += matrix[j];
            }
        }
        printMatrix(size, sumMatrix);
    }

    public static void checkAlgorithm() {
        final int size = 8; 
        final int count = 2215764;
        final int divisor = 122;
        final int sumMatrix[] = new int[size * size];
        for (int i = 0; i < count/divisor ; i++) {
            final int[] matrix = generate(size);
            for (int j = 0; j < sumMatrix.length; j++) {
                sumMatrix[j] += matrix[j];
            }
        }
        int total = 0;
        for(int i=0; i < sumMatrix.length; i++) {
            total += sumMatrix[i];
        }
        final double factor = (double)total / (count/divisor);
        System.out.println("Factor=" + factor + " (theory=16.0)");
    }

    public static void benchmark(final int size, final int count,
            final boolean parallel) {
        final long begin = System.currentTimeMillis();
        if (!parallel) {
            for (int i = 0; i < count; i++) {
                generate(size);
            }
        } else {
            IntStream.range(0, count).parallel().forEach(i -> generate(size));
        }
        final long end = System.currentTimeMillis();
        System.out.println("rate="
                + (double) (end - begin) / count + "ms/matrix");
    }

    public static void main(String[] args) {
        checkAlgorithm();
        benchmark(8, 10000, true);
        //printStatistics(8, 2215764/36);
        printStatistics(8, 2215764);
    }
}

Вывод:

Factor=16.0 (theory=16.0)
rate=0.2835ms/matrix
552969 554643 552895 554632 555680 552753 554567 553389 
554071 554847 553441 553315 553425 553883 554485 554061 
554272 552633 555130 553699 553604 554298 553864 554028 
554118 554299 553565 552986 553786 554473 553530 554771 
554474 553604 554473 554231 553617 553556 553581 553992 
554960 554572 552861 552732 553782 554039 553921 554661 
553578 553253 555721 554235 554107 553676 553776 553182 
553086 553677 553442 555698 553527 554850 553804 553444

Ответ 5

Вот очень быстрый подход к созданию матрицы row by row, написанной на Java:

public static void main(String[] args) throws Exception {
    int n = 100;
    Random rnd = new Random();
    byte[] mat = new byte[n*n];
    byte[] colCount = new byte[n];

    //generate row by row
    for (int x = 0; x < n; x++) {               

        //generate a random first bit
        int b1 = rnd.nextInt(n);
        while ( (x > 0 && mat[(x-1)*n + b1] == 1) ||    //not adjacent to the one above
                (colCount[b1] == 2)                     //not in a column which has 2
                ) b1 = rnd.nextInt(n);


        //generate a second bit, not equal to the first one
        int b2 = rnd.nextInt(n);
        while ( (b2 == b1) ||                           //not the same as bit 1
                (x > 0 && mat[(x-1)*n + b2] == 1) ||    //not adjacent to the one above
                (colCount[b2] == 2) ||                  //not in a column which has 2
                (b2 == b1 - 1) ||                       //not adjacent to b1
                (b2 == b1 + 1)
                ) b2 = rnd.nextInt(n);


        //fill the matrix values and increment column counts
        mat[x*n + b1] = 1;
        mat[x*n + b2] = 1;              
        colCount[b1]++;
        colCount[b2]++;
    }

    String arr = Arrays.toString(mat).substring(1, n*n*3 - 1);
    System.out.println(arr.replaceAll("(.{" + n*3 + "})", "$1\n"));     
}

Он по существу генерирует каждую случайную строку за раз. Если строка будет нарушать какое-либо из условий, она будет сгенерирована снова (опять же случайным образом). Я считаю, что это также удовлетворит условие 4.

Добавление быстрой заметки о том, что она будет вращаться навсегда для N-s, где нет решений (например, N = 3).