Определяющий коэффициент члена x^m в (x^2 + x + 1)^n четен или нечетен

Для заданных целых чисел n и m определить, является ли коэффициент x^m члена в (x^2+x+1)^n четным или нечетным?

Например, если n=3 и m=4, (x^2+x+1)^3 = x^6 + 3x^5 + [[6x^4]] + 7x^3 + 6x^2 + 3x + 1, значит, коэффициент x^4 члена равен 6 (=четное).

n и m равно 10^12, и я хочу вычислить за несколько секунд, поэтому вы не можете вычислить за линейное время.

У вас есть какой-нибудь эффективный алгоритм?


person square1001    schedule 29.04.2017    source источник
comment
Что вы имеете в виду под коэффициентом x^m?   -  person trincot    schedule 29.04.2017
comment
@trincot Добавлен пример.   -  person square1001    schedule 29.04.2017
comment
@VincentvanderWeele Правда? Это не все математика. Это вопрос программирования, потому что я хочу рассчитать БЫСТРЕЕ, а не 10 или 20 минут в одной паре значений. Это вопрос алгоритма.   -  person square1001    schedule 29.04.2017
comment
Если вопрос задан таким образом, вы можете быть уверены, что вам не нужно вычислять фактический коэффициент. Я бы, вероятно, начал с преобразования всего в двоичный формат и поискал, какие термины я могу сразу исключить.   -  person biziclop    schedule 29.04.2017


Ответы (4)


Да, линейное время по количеству битов на входе.

Рассматриваемые коэффициенты являются трехчленными коэффициентами T(n, m). Для биномиальных коэффициентов мы будем использовать теорему Лукаса; давайте разработаем трехчленный аналог для p = 2.

Работая mod 2 и следуя доказательству Натана Файна,

(1 + x + x^2)^{2^i} = 1 + x^{2^i} + x^{2^{2 i}}

(1 + x + x^2)^n
    = prod_i ((1 + x + x^2)^{2^i n_i})
        where n = sum_i n_i 2^i and n_i in {0, 1} for all i
        (i.e., n_i is the binary representation of n
    = prod_i (1 + x^{2^i n_i} + x^{2^i 2 n_i})
    = prod_i sum_{m_i = 0}^{2 n_i} x^{2^i}
    = sum_{(m_i)} prod_i x^{2^i m_i}
        taken over sequences (m_i) where 0 ≤ m_i ≤ 2 n_i.

В биномиальном случае следующим шагом является наблюдение, что для коэффициента x^m существует не более одного выбора (m_i), чьи x^{2^i m_i} множители имеют правильное произведение, т. е. двоичное представление m.

В трехчленном случае мы должны рассматривать двоичные псевдопредставления (m_i) из m, где псевдобиты могут быть равны нулю, единице или двум. Вклад в сумму есть тогда и только тогда, когда для всех i таких, что n_i = 0, имеем m_i = 0.

Мы можем написать автомат, который сканирует n и m бит за битом. Состояние a является начальным и принимающим.

a (0:0:nm') -> a nm'    [emit 0]
a (1:0:nm') -> a nm'    [emit 0]
            -> b nm'    [emit 2]
a (1:1:nm') -> a nm'    [emit 1]

b (0:1:nm') -> a nm'    [emit 0]
b (1:0:nm') -> b nm'    [emit 1]
b (1:1:nm') -> a nm'    [emit 0]
            -> b nm'    [emit 2]

Мы можем использовать динамическое программирование для подсчета путей. В виде кода:

def trinomial_mod_two(n, m):
    a, b = 1, 0
    while m:
        n1, n = n & 1, n >> 1
        m1, m = m & 1, m >> 1
        if n1:
            if m1:
                a, b = a ^ b, b
            else:
                a, b = a, a ^ b
        elif m1:
            a, b = b, 0
        else:
            a, b = a, 0
    return a

Версия без веток для хихиканья:

def trinomial_mod_two_branchless(n, m):
    a, b = 1, 0
    while m:
        n1, n = n & 1, n >> 1
        m1, m = m & 1, m >> 1
        a, b = ((n1 | ~m1) & a) ^ (m1 & b), ((n1 & ~m1) & a) ^ (n1 & b)
    return a
person David Eisenstat    schedule 29.04.2017
comment
Очень умно, и я хотел бы подумать об этом подходе. Я думаю, что должен быть более ясный способ объяснить это, хотя. Кое-что о подсчете способов построить m из сложения некоторой комбинации битов n (или тех битов, которые были сдвинуты влево один раз). - person Paul Hankin; 29.04.2017
comment
Вау! Очень красивое решение! - person square1001; 30.04.2017

Во-первых, заметим, что если кого-то интересует только, является ли коэффициент при x^m нечетным или четным, то можно рассматривать коэффициенты многочлена как элементы конечного поля F2.

Обратите внимание на (1+x+x^2)^2 = (1+x^2+x^4) mod 2, потому что все перекрестные члены четные. На самом деле, если n является степенью числа 2, то (1+x+x^2)^n = (1 + x^n + x^2n) mod 2.

Для общего n запишите его как сумму степеней двойки (то есть в двоичном формате).

n = (2^a1 + 2^a2 + 2^a3 + ... + 2^ak).

Затем перемножьте степени, соответствующие каждой степени числа 2:

(1+x+x^2)^n = (1+x^(2^a1)+x^(2^(a1+1))) * ((1+x^(2^a2)+x^(2^(a2+1))) * ...

Каждое из условий в этом произведении теперь имеет только 3 множителя, и их максимум 35 или 36, если n ограничено 10 ^ 12. Так что их легко умножить вместе.

Вот некоторый код Python, который делает это:

# poly_times computes the product of polynomials
# p and q over the field F2. They are each
# represented by a set of non-zero coefficients.
# For example set([0, 2, 5]) corresponds to x^0 + x^2 + x^5.
def poly_times(p, q):
    result = set()
    for i in p:
        for j in q:
            if i+j in result:
                result.remove(i+j)
            else:
                result.add(i+j)
    return result

# Return the coefficient of x^m in (1+x+x^2)^n.
def coeff(n, m):
    prod = set([0])
    i = 0
    while n:
        if n % 2:
            prod = poly_times(prod, [0, 2**i, 2**(i+1)])
        i += 1
        n >>= 1
    return 1 if m in prod else 0

for m in xrange(10):
    print m, coeff(3, m)

print coeff(1947248745, 1947248745034)    

Оптимизация

Для n с большим количеством установленных битов это становится слишком медленным, когда n приближается к 10^12.

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

# poly_times computes the product of polynomials
# p and q over the field F2. They are each
# represented by a set of non-zero coefficients.
# For example set([0, 2, 5]) corresponds to x^0 + x^2 + x^5.
# Coefficients larger than m are discarded.
def poly_times(p, q, m):
    result = set()
    for i in p:
        for j in q:
            if i + j > m:
                continue
            if i+j in result:
                result.remove(i+j)
            else:
                result.add(i+j)
    return result

# Return the coefficient of x^m in (1+x+x^2)^n.
def coeff(n, m):
    if m > 2*n: return 0
    prod = [set([0]), set([0])]
    i = 0
    while n:
        if n % 2:
            prod[i//20] = poly_times(prod[i//20], [0, 2**i, 2**(i+1)], m)
        i += 1
        n >>= 1
    s = 0
    for x in prod[0]:
        s += m-x in prod[1]
    return s % 2

for m in xrange(10):
    print m, coeff(3, m)

print coeff(0xffffffffff, 0xffffffffff)

Обратите внимание, что это может вычислить coeff(0xffffffffff, 0xffffffffff) за несколько секунд, а 0xffffffffff больше, чем 10**12.

person Paul Hankin    schedule 29.04.2017
comment
Хак кейс: n=2^35-1, m=2^35-1 (думаю количество нечетных коэффициентов будет много) - person square1001; 29.04.2017
comment
В худшем случае это определенно линейно, но его довольно сложно вычислить, поскольку, если n имеет много установленных битов, то есть много общих терминов. Но я думал, что вопрос в том, чтобы найти решение, работающее за несколько секунд, а не в вычислительной сложности? Я догадался из вопроса, что это был вопрос типа онлайн-судьи. Если это недостаточно быстро, есть много возможностей для дальнейшей оптимизации этого решения, даже помимо пары предложений, которые я сделал в вопросе. - person Paul Hankin; 29.04.2017
comment
Я думаю, вы правы, что это слишком медленно, когда n имеет много битов. - person Paul Hankin; 29.04.2017
comment
У меня есть оптимизация, которая делает его очень быстрым. - person Paul Hankin; 29.04.2017
comment
Хорошо, теперь для n=2^40-1, m=2^40-1 требуется 4 секунды. - person Paul Hankin; 29.04.2017

Процентный коэффициент зависит от количества способов, которыми можно выбрать n терминов из x² + x + 1 так, чтобы сумма степеней выбранных терминов была равна м. Эти способы могут быть сгруппированы в группы с одинаковым количеством выбранных терминов и терминов x (количество выбранных раз 1 следует из тот).

Пусть a будет количеством терминов , b количеством терминов x, а c количество терминов 1 в определенной группе.

Тогда выполняются следующие равенства:

2a + b = m
a + b + c = n

Очевидно, что в целом существует несколько групп с разными значениями для a, b, c. После определения a также определяются значения для b и c. Таким образом, нужно только перебрать возможные значения для a, чтобы получить все группы.

Если бы вы написали алгоритм грубой силы для получения самого коэффициента, на Python он выглядел бы так:

def combi(n, k): # number of ways to take k elements from n elements
    import math
    f = math.factorial
    return f(n) // f(k) // f(n-k)    

def get_coeff(n, m):
    if m > n * 2 or n < 0 or m < 0: # basic argument check
        return None
    if m > n: # mirrored situation is the same
        m = 2*n - m            
    coeff = 0
    for a in range(0, m//2+1):
        b = m - 2*a
        coeff += combi(n, a) * combi(n-a, b)
    return coeff

Функция combi(n, k) вернет количество способов взять k элементов из n элементов, т. е. биномиальный коэффициент.

Произведение двух вызовов combi отвечает на следующий вопрос:

Сколькими способами я могу взять a, умноженное на , и b, умноженное на x? Обратите внимание, что количество способов, которыми может быть выбран постоянный член, равно 1 после того, как были сделаны 2 других выбора.

Теперь, поскольку нам нужно только знать, является ли конечный коэффициент нечетным или четным, нам также нужно только знать, является ли биномиальный коэффициент нечетным или четным. Как объясняется на math.stackexchange.com, это можно определить с помощью простого бита операция:

def combi_parity(n, k):
    return (k & (n-k)) == 0

def get_coeff_parity(n, m):
    if m > n * 2 or n < 0 or m < 0: # basic argument check
        return None
    if m > n:
        m = 2*n - m # mirrored situation is the same
    coeff_odd = 0
    for a in range(0, m//2+1):
        b = m - 2*a
        # A product is odd only when both factors are odd (so we perform an "and")
        # A sum changes parity whenever the added term is odd (so we perform a "xor")
        coeff_odd ^= combi_parity(n, a) and combi_parity(n-a, b) 
    return coeff_odd

Посмотрите, как он работает на repl.it.

person trincot    schedule 29.04.2017
comment
Проголосовал. Я тоже попробовал этот подход, но не смог сделать это достаточно быстро для больших m, n. Я даже немного оптимизировал, чтобы a колебался только в числах с combi_parity(n, a)==1. Вот мой код C, который все еще работает на моем ноутбуке: gist.github.com/paulhankin/5f5aa05231ce8659b6b040d849fd800e . - person Paul Hankin; 29.04.2017
comment
Ваша формулировка c(n, k) эквивалентна, но лучше той, которую использовал я: n&k==k. Я собираюсь украсть его для своей программы! - person Paul Hankin; 29.04.2017
comment
Вы не используете большие числа в своем тестовом примере: 2 ^ 40-1 должно быть 2 ** 40-1. Извините, я перепутал питон и математическую нотацию. 2^40-1 равно 37;( - person Paul Hankin; 29.04.2017
comment
Ааааааааа, спасибо! Я обновил образец вызова в указанной ссылке. - person trincot; 29.04.2017
comment
Мой вложенный цикл на самом деле представляет собой (попытку) оптимизации, которая выполняет a (который в моем коде называется k) только для чисел, для которых combi_parity(n, a) имеет значение True. - person Paul Hankin; 29.04.2017
comment
Хорошо, я пропустил это. - person trincot; 29.04.2017

Хорошо, я просто придумал решение. Вот оно:

  1. Представьте, что уравнение записано n раз, (a.x^2 + b.x^1 + c).(a.x^2 + b.x^1 + c)...n раз. a, b и c - константы, которые я принял вообще.
  2. Теперь нам нужно выбрать из каждого термина так, чтобы результат умножения всех таких терминов давал x^m
  3. Теперь я могу сказать, что нам нужно найти решения уравнения t1.2 + t2 = m, где t1 - это отсутствие появления x^2 и t2 из x. Это связано с тем, что t1 и t2 в этом случае будут иметь форму k.x^m (k является константой). Это нахождение интегральных диофантовых решений этого уравнения, то есть нахождение всех удовлетворяющих пар {t1, t2}
  4. Теперь нам нужно применить небольшую перестановку для нахождения здесь коэффициента. Предполагая, что у вас есть одно из решений {1, 2} для предыдущего шага, тогда для этой диады коэффициент будет (1^1.nC1).(2^2.(n-1)C2), который будет одним из составляющих коэффициента. Если вы суммируете все такие коэффициенты, соответствующие всем диофантовым решениям, вы получите коэффициент.

Реализация вышеизложенного алгоритмически заняла бы у меня некоторое время, поэтому я опубликовал шаги.

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

РЕДАКТИРОВАТЬ: Как было задано для примера,

  1. Допустим, у нас есть уравнение (x^2 + x^1 + x^1)^3, и мы должны найти коэффициент x^3. Таким образом, у нас есть m = 3.
  2. Отдельно писать уравнение для того, чтобы визуально видеть шаг. Это,

    (x^2 + x^1 + x^1).(x^2 + x^1 + x^1).(x^2 + x^1 + x^1)

  3. Теперь мы хотим выбрать любой из {x^2, x^1, 1} из каждого из них. Будет несколько способов выбрать его, чтобы сделать умножение формы, x^3

  4. Чтобы решить эту проблему, мы можем написать уравнение 2.a + b = 3, где a равно количеству раз, выбранных x^2, а b равно количеству раз, когда x выбрано. Решения этого уравнения: {0, 3} и {1, 1}. Теперь, поскольку мы также должны учитывать порядок, в котором мы их выбираем, мы применим комбинаторику на следующем шаге.

  5. Коэффициент будет 2^0.3C0.3^3.3C3 + 2^1.3C1.3^1.2C1. Теперь здесь, в первом члене, 2^0.3C0.3^3.3C3, 3C0 означает выбор 0 вхождений x^2 с 3C3 означает 3 вхождения x, что даст x^3, но мы также умножаем на 2^0, потому что 2 является коэффициентом x^2 в уравнении, и аналогично, 3^3 потому что 3 является коэффициентом x. Точно так же вы подходите ко второму члену, соответствующему {1, 1}.

  6. В сумме получается 63, что вы можете проверить, умножив вручную, и вы получите 63.

Надеюсь, я ясно.

person Manish Kumar Sharma    schedule 29.04.2017
comment
Я не понимаю. Можете ли вы привести пример? - person square1001; 29.04.2017
comment
@square1001: я готовлю пример. Дай мне пару минут. - person Manish Kumar Sharma; 29.04.2017
comment
@square1001: см. пример. - person Manish Kumar Sharma; 29.04.2017
comment
Но если n и m велики, вы не можете использовать грубую силу {a, b}. - person square1001; 29.04.2017
comment
Подводя итог этому методу: вычислить sum(c(n, k)c(n-k, m-2k) for k=0..m/2). Вы можете сделать это довольно быстро, вычислив его по модулю 2 и отметив, что c(a, b)=1 по модулю 2, только если a&b==b (так что вам нужно только перебирать k с битами подмножества битов n) . Но мне не удалось заставить C-версию работать достаточно быстро, когда n и m большие. Если вы хотите поиграть с моим кодом: gist.github.com/paulhankin/5f5aa05231ce8659b6b040d849fd800e - person Paul Hankin; 29.04.2017