Подтвердить что ты не робот

Быстрый способ найти нулевой ответ

Для каждого массива длины n + h-1 со значениями от 0 и 1, я хотел бы проверить, существует ли другой ненулевой массив длины n со значениями от -1,0,1, так что все h внутренние произведения равны нулю. Мой наивный способ сделать это -

import numpy as np
import itertools
(n,h)= 4,3
for longtuple in itertools.product([0,1], repeat = n+h-1):
    bad = 0
    for v in itertools.product([-1,0,1], repeat = n):
        if not any(v):
            continue
        if (not np.correlate(v, longtuple, 'valid').any()):
            bad = 1
            break
    if (bad == 0):
        print "Good"
        print longtuple

Это очень медленно, если мы установили n = 19 и h = 10, что я хотел бы проверить.

Моя цель - найти единственный "хороший" массив длины n+h-1. Есть ли способ ускорить это так, чтобы выполнялись n = 19 и h = 10?

Нынешний наивный подход принимает 2 ^ (n + h-1) 3 ^ (n) итераций, каждый из которых занимает примерно n времени. Это 311,992,186,885,373,952 итераций для n = 19 и h = 10, что невозможно.

Примечание 1 Изменено convolve до correlate, так что код считает v правильным путем.


10 июля 2015

Проблема по-прежнему открыта без решения, достаточно быстрого для n=19 и h=10.

4b9b3361

Ответ 1

Рассмотрим следующий подход "встретить в середине".

Во-первых, переформулируйте ситуацию в матричной формулировке, предоставляемой leekaiinthesky.

Далее отметим, что нам нужно рассматривать только "короткие" векторы s формы {0,1}^n (т.е. короткие векторы, содержащие только 0 и 1), если мы изменим задачу на поиск h x n Hankel matrix H от 0 и 1, что Hs1 никогда не будет равно Hs2 для двух разных коротких векторов 0 и 1. Это связано с тем, что Hs1 = Hs2 означает H(s1-s2)=0, что означает, что существует вектор v 1, 0 и -1, а именно s1-s2, такой, что Hv = 0; наоборот, если Hv = 0 для v в {-1,0,1}^n, то мы можем найти s1 и s2 в {0,1}^n такие, что v = s1 - s2 so Hs1 = Hs2.

Когда n=19 есть только 524 288 векторов s в {0,1}^n, чтобы попробовать; хэш результатов Hs, и если один и тот же результат произойдет дважды, то H не подходит и попробуйте другой H. С точки зрения памяти этот подход вполне возможен. Существуют 2^(n+h-1) ганкелевые матрицы H; когда n=19 и h=10 - 268 435 456 матриц. Тест 2^38, или 274 877 906 944, каждый из которых имеет около nh операций, чтобы умножить матрицу H и вектор s, около 52 триллионов операций. Это кажется выполнимым, нет?

Поскольку теперь вы используете только 0 и 1, а не -1, вы также можете ускорить процесс, используя бит-операции (shift, and, and count 1's).

Обновление

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

Я получил 00000000010010111000100100 для n = 17 и h = 10 через несколько минут, а 000000111011110001001101011 для n = 18 и h = 10 через некоторое время. Я собираюсь запустить его для n = 19 и h = 10.

#include <iostream>
#include <bitset>
#include <unordered_set>

/* Count the number of 1 bits in 32 bit int x in 21 instructions.
 * From /Hackers Delight/ by Henry S. Warren, Jr., 5-2
 */
int count1Bits(int x) {
  x = x - ((x >> 1) & 0x55555555);
  x = (x & 0x33333333) + ((x >> 2) & 0x33333333);
  x = (x + (x >> 4)) & 0x0F0F0F0F;
  x = x + (x >> 8);
  x = x + (x >> 16);
  return x & 0x0000003F;
}

int main () {
  const int n = 19;
  const int h = 10;
  std::unordered_set<long> dotProductSet;

  // look at all 2^(n+h-1) possibilities for longVec
  // upper n bits cannot be all 0 so we can start with 1 in pos h
  for (int longVec = (1 << (h-1)); longVec < (1 << (n+h-1)); ++longVec) {

    dotProductSet.clear();
    bool good = true;

    // now look at all n digit non-zero shortVecs
    for (int shortVec = 1; shortVec < (1 << n); ++shortVec) {

      // longVec dot products with shifted shortVecs generates h results
      // each between 0 and n inclusive, can encode as h digit number in
      // base n+1, up to (n+1)^h = 20^10 approx 13 digits, need long
      long dotProduct = 0;

      // encode h dot products of shifted shortVec with longVec
      // as base n+1 integer
      for(int startShort = 0; startShort < h; ++startShort) {
        int shortVecShifted = shortVec << startShort;
        dotProduct *= n+1;
        dotProduct += count1Bits(longVec & shortVecShifted);
      }

      auto ret = dotProductSet.insert(dotProduct);
      if (!ret.second) {
        good = false;
        break;
      }
    }

    if (good) {
      std::cout << std::bitset<(n+h-1)>(longVec) << std::endl;
      break;
    }
  }

  return 0;
}

Второе обновление

Программа для n = 19 и h = 10 работала в течение двух недель в фоновом режиме на моем ноутбуке. В конце он просто вышел без печати каких-либо результатов. Запрет какой-либо ошибки в программе, похоже, что нет длинных векторов с требуемым свойством. Я предлагаю искать теоретические причины, почему таких длинных векторов нет. Возможно, какой-то счетный аргумент будет работать.

Ответ 2

Может быть более быстрый способ *

То, что вы ищете, связано с концепцией ядро ​​ или пустое пространство матрицы.

В частности, для каждого n + h-1 "longtuple" и заданного n построим h по n матрице, строки которой являются n-подмножествами длинной вершины. Другими словами, если ваш longtuple [0,0,0,1,0,0] и n = 3, то ваша матрица:

[[0 0 0]
 [0 0 1]
 [0 1 0]
 [1 0 0]]

Вызвать эту матрицу A. Вы ищете вектор x, для которого Ax= 0, где 0 - вектор всех 0s, Если такой x существует (то есть не все все 0s) и может быть масштабирован, чтобы содержать только {-1, 0, 1}, то вы хотите выбросить A и перейдите к следующему длинному набору.

Я не совсем уверен, что (теоретически наиболее эффективная) вычислительная сложность вычисления ядра, но он, кажется, имеет порядок O (h + n) ^ 3 или так, что в любом случае a намного лучше, чем O (3 ^ n). См. Ссылку Wikipedia выше или Python (NumPy, SciPy), нахождение нулевого пространства матрицы для некоторых примеров того, как вычислить ядро.

В любом случае, как только вы идентифицируете ядро, вам нужно будет сделать дополнительную работу, чтобы выяснить, существуют ли какие-либо векторы формы {-1, 0, 1} ^ n, но я не думаю, что как большая вычислительная нагрузка.

* NB. В комментариях @vib указывает, что на самом деле это может быть большой вычислительной нагрузкой. Я не уверен, какой лучший алгоритм для выяснения, пересекаются ли эти векторы ядра. Возможно, это невозможно решить за полиномиальное время, и в этом случае этот ответ не обеспечивает ускорения исходной проблемы!

Пример кода

Адаптация кода из другого вопроса о переполнении стека, связанного выше, для примера, который вы указали в комментариях:

import scipy
from scipy import linalg, matrix
def null(A, eps=1e-15):
    u, s, vh = scipy.linalg.svd(A)
    null_mask = (s <= eps)
    null_space = scipy.compress(null_mask, vh, axis=0)
    return scipy.transpose(null_space)

A = matrix([[0,0,0,1],[0,0,1,0],[0,1,0,0],[0,0,0,0]])
print null(A)

#> [[-1.]
#>  [ 0.]
#>  [ 0.]
#>  [ 0.]]

В коде приведен пример (фактически тот же пример, который вы дали) n-кортежа, который делает недействительным [0, 0, 0, 1, 0, 0] как "хорошую" длинную строку. Если код возвратил [], то, предположительно, нет такого n-кортежа, а longtuple "хорош". (Если код действительно что-то возвращает, вам все равно нужно проверить часть {-1, 0, 1}.)

Дальнейшие мысли

Существует ли такой x, без учета ограничения {-1, 0, 1}, эквивалентен вопросу о том, является ли недействительность (размерность ядро) A больше 0. Это было бы эквивалентно заданию того, равен ли ранг A n. Поэтому, если вы нашли какой-то способ быть умным в отношении ограничения {-1, 0, 1} и сломать его просто для того, чтобы вычислить ранг A, я уверен, что это можно было бы сделать еще быстрее.

Кстати, очень вероятно, что вы (или тот, кто дал вам эту проблему), могут знать все это уже... В противном случае, почему бы вы назвали длину длинной строки "n + h-1", если вы еще не началось с матрицы высоты h...!

Ответ 3

Это лишь частичный ответ, так как это все еще кажется слишком медленным, чтобы проверить случай n=19, h=10 (и в этом случае нет даже хороших векторов).

Вот реализация алгоритма проверки, который @vib описывает и использует произвольную выборку по всем векторам 2^(n+h-1), в Mathematica.

TestFullNullSpace[A_, ns_] := Module[{dim, n},
  {dim, n} = Dimensions[ns];
  For[i = 1, i < 3^dim, i++,
   testvec = Mod[IntegerDigits[i, 3, dim].ns, 3] /. {2 -> -1};
   If[Norm[A.testvec] == 0,
    Return[False]];
   ];
  Return[True];
]

n = 17;
h = 10;

Do[
 v = Table[RandomChoice[{0, 1}], {n + h - 1}];
 A = Table[v[[i ;; i + n - 1]], {i, 1, h}];
 ns = NullSpace[A, Modulus -> 3] /. {2 -> -1};
 If[TestFullNullSpace[A, ns],
  Print[v]];,
 {1000}]

Пример вывода для вышеуказанного запуска через несколько секунд вычисления:

{0,0,1,1,0,0,0,0,0,0,1,0,1,1,1,0,1,0,1,1,0,0,0,1,1,0}
{1,1,0,1,0,0,0,1,1,0,1,1,1,1,1,0,1,0,1,0,0,1,1,0,0,0}
{1,0,1,1,1,1,1,0,0,0,1,1,0,1,0,1,1,0,0,1,1,0,1,1,1,0}
{0,0,0,0,1,0,1,1,1,0,1,1,0,0,1,1,1,1,0,1,0,1,0,0,1,1}

Итак, из 1000 проверенных векторов 4 были "хорошими" (если у меня нет ошибки). К сожалению, для n=18 я провел это несколько минут и до сих пор не нашел "хорошего" вектора. Я не знаю, существуют ли они или просто редко встречаются.

Ответ 4

Ниже приведен алгоритм, который уменьшает сложность от 3 ^ n до 3 ^ {n-h}.

Пусть v_1, v_2,.., v_h - векторы, которые вам нужно ортогонально.

Рассмотрим векторное пространство (Z/3Z) ^ n. Пусть v'_1,..., v'_h - естественные включения v_1,..., v_h в этом пространстве.

Пусть теперь w - вектор с коэффициентами в {-1,0,1}, w '- вектор (Z/3Z) ^ n, полученный естественным образом, рассматривая w как вектор (Z/3Z) ^ п. Тогда необходимым условием для w иметь нулевое скалярное произведение с v_1,..., v_h (в R) является то, что w 'имеет нулевое скалярное произведение (в (Z/3Z) ^ n) с v'_1,..., v "_h.

Теперь вы можете легко определить w ', которые имеют нулевое скалярное произведение с v'_1,..., v'_h. Они образуют пространство размером 3 ^ {n-h}. Затем вам нужно проверить для каждого из них, был ли связанный w фактически ортогонален всем v_i.

Ответ 5

Вот такой подход, который сводит его к O(n*h*3^(n/2 + 1)). Это плохо масштабируется, но достаточно хорошо для вашего случая использования.

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

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

Не забывайте игнорировать ответ, который является всем 0s!