Проверить быстрее, чем вычислить: Фрейвалдс и тождество многочленов

Если ты запишешь проверку Фрейвалдса, объяснишь, откуда в ней 1/2, и скажешь, почему проверять многочлены надо в случайной точке и по простому модулю, - листай до «Теперь сам».

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

Сторонний сервис перемножил две матрицы 1000 × 1000 и вернул результат. Хочется убедиться, что он не ошибся. Перемножить самому - миллиард умножений.

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

На что это похоже

Сверка двух длинных ведомостей. Сравнивать их строка за строкой долго. Вместо этого обе стороны складывают свои суммы, но не все подряд, а только по строкам, которые выпали по жребию - один и тот же жребий для обеих. Если ведомости одинаковы, суммы совпадут при любом жребии. Если отличаются хотя бы в одной строке, то эта строка с вероятностью 1/2 попала в жребий - и тогда суммы, скорее всего, разойдутся.

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

Механизм: случайный вектор и случайная точка

Проверка Фрейвалдса. Выбираем вектор r длины n, каждая координата - 0 или 1 с вероятностью 1/2 независимо. Считаем A · (B · r) и C · r и сравниваем.

B · r        матрица на вектор:   n² действий, получается вектор
A · (B · r)  матрица на вектор:   n² действий
C · r        матрица на вектор:   n² действий

Скобки здесь - вся хитрость: A · (B · r) равно (A · B) · r, но считается без произведения матриц.

Если A · B = C, векторы совпадут при любом r. Пусть A · B ≠ C. Обозначим D = A · B − C; это ненулевая матрица, и проверка ошибается, когда D · r всё равно оказался нулевым вектором. В D есть ненулевой элемент; пусть он стоит в строке i и столбце j. Координата i вектора D · r:

dᵢ₁·r₁ + … + dᵢⱼ·rⱼ + … + dᵢₙ·rₙ          где dᵢⱼ ≠ 0

Зафиксируем все координаты r, кроме rⱼ. Остальные слагаемые дают какое-то число s, и вся сумма равна s + dᵢⱼ · rⱼ. При rⱼ = 0 она равна s, при rⱼ = 1 - s + dᵢⱼ. Эти два числа различны, значит нулём может быть не более чем одно из них. А rⱼ выбирается монеткой независимо от остальных. Поэтому с вероятностью не меньше 1/2 координата i ненулевая, и ошибка поймана.

Тождество многочленов. Даны два выражения, задающие многочлены, например (x − 1)(x − 2)(x − 3)(x − 4) и x⁴ − 10x³ + 35x² − 50x + 24. Равны ли они как многочлены? Раскрывать скобки - d² действий для степени d. Вместо этого берём случайное x и вычисляем оба выражения.

Опора - факт из алгебры: ненулевой многочлен степени d имеет не больше d корней. Разность двух многочленов степени не выше d, разных по модулю p, - ненулевой многочлен степени не выше d, и значения исходных совпадают только в его корнях. Если x выбрано из p возможных значений, вероятность попасть в корень не больше d / p.

Чтобы числа не росли, всё считают по модулю простого числа p, а x берут от 0 до p − 1. Простота модуля существенна: факт о числе корней верен там, где на любое ненулевое число можно делить, а по составному модулю это не так.

Правило

Проверка Фрейвалдса: для случайного вектора r из нулей и единиц сравнить A · (B · r) с C · r. Цена - 3n². Если A · B = C, ответ всегда «равно»; если нет, ответ «равно» выдаётся с вероятностью не больше 1/2 за одну проверку и не больше 1/2ᵏ за k независимых. Проверка тождества многочленов степени не выше d в случайной точке по простому модулю p ошибается с вероятностью не больше d / p. Обе проверки - алгоритмы Монте-Карло с односторонней ошибкой: ответу «не равно» можно верить всегда.

Разбор: сорок проверок верного и испорченного произведения

Матрицы 12 × 12 со случайными элементами от −9 до 9. c - настоящее произведение, посчитанное тремя вложенными циклами. Сначала сорок проверок на верном c, затем один элемент c увеличен на единицу и сделано ещё сорок.

import random

def umnozhit_na_vektor(m, v):
    n = len(v)
    otvet = []
    for i in range(n):
        s = 0
        for j in range(n):
            s = s + m[i][j] * v[j]
        otvet.append(s)
    return otvet

def proverka(a, b, c):
    n = len(a)
    r = []
    for i in range(n):
        r.append(random.randint(0, 1))
    return umnozhit_na_vektor(a, umnozhit_na_vektor(b, r)) == umnozhit_na_vektor(c, r)

random.seed(34)
n = 12
a = []
b = []
for i in range(n):
    stroka_a = []
    stroka_b = []
    for j in range(n):
        stroka_a.append(random.randint(-9, 9))
        stroka_b.append(random.randint(-9, 9))
    a.append(stroka_a)
    b.append(stroka_b)

c = []
for i in range(n):
    stroka = []
    for j in range(n):
        s = 0
        for k in range(n):
            s = s + a[i][k] * b[k][j]
        stroka.append(s)
    c.append(stroka)

vernyh = 0
for opyt in range(40):
    if proverka(a, b, c):
        vernyh = vernyh + 1

c[3][7] = c[3][7] + 1
propuskov = 0
for opyt in range(40):
    if proverka(a, b, c):
        propuskov = propuskov + 1
print(vernyh, propuskov)

40 22

Верное произведение прошло все сорок проверок - иначе и быть не может. Испорченное прошло 22 из сорока: ошибка в одном элементе c[3][7] видна, только когда r[7] = 1, а это половина случаев. Одна проверка ненадёжна, но вероятность того, что испорченная матрица пройдёт двадцать проверок подряд, - 1/2²⁰, одна миллионная.

Повторов можно сделать меньше, если брать координаты r не из двух значений, а из тысячи: то же рассуждение даёт, что сумма s + dᵢⱼ · rⱼ обращается в ноль не более чем при одном значении rⱼ из тысячи, и одна проверка ошибается с вероятностью не больше 0,001.

Разбор: случайная точка и точки, которые выбрал бы человек

Многочлен (x − 1)(x − 2)(x − 3)(x − 4) и два набора коэффициентов: верный и с ошибкой в двух местах. proizvedenie перемножает скобки, gorner вычисляет многочлен по коэффициентам схемой Горнера - коэффициенты идут от старшей степени к свободному члену; всё по модулю простого числа 67 108 859.

import random

p = 67108859

def proizvedenie(korni, x):
    y = 1
    for r in korni:
        y = y * (x - r) % p
    return y

def gorner(koef, x):
    y = 0
    for k in koef:
        y = (y * x + k) % p
    return y

korni = [1, 2, 3, 4]
verno = [1, -10, 35, -50, 24]
oshibka = [1, -9, 35, -51, 24]

random.seed(35)
x = random.randrange(p)
print(proizvedenie(korni, x) == gorner(verno, x), proizvedenie(korni, x) == gorner(oshibka, x))
for x in [0, 1, p - 1]:
    print(x, proizvedenie(korni, x) == gorner(oshibka, x))

True False
0 True
1 True
67108858 True

Случайная точка подтвердила верные коэффициенты и отвергла ошибочные. А три точки, которые первыми приходят в голову, - ноль, единица и минус единица (по модулю это p − 1) - ошибку пропустили все три. Ошибочный многочлен отличается от верного на x³ − x, а у этой разности корни как раз 0, 1 и −1.

Это и есть причина, по которой точка должна быть случайной. Разность двух многочленов степени четыре имеет не больше четырёх корней из 67 миллионов возможных точек. Проверяющий, который всегда подставляет одни и те же «удобные» числа, гарантированно пропускает ошибки, чья разность зануляется именно в них; против случайной точки такую ошибку подобрать нельзя. Вероятность промаха здесь не больше 4 / 67 108 859.

Граница приёма - модуль. По составному модулю корней бывает больше степени: у x² − 1 по модулю 8 их четыре - 1, 3, 5 и 7. Оценка d / p опирается на простоту p.

Что ломается без этого

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

Дробные числа. Для матриц из чисел с плавающей точкой A · (B · r) и C · r различаются погрешностями округления даже при верном C. Сравнивать нужно с допуском, и рассуждение про «ровно один корень» превращается в приближённое.

Ответу «равно» верят как доказательству. Проверка односторонняя: «не равно» - это найденное расхождение, а «равно» - отсутствие расхождения в k попытках.

Переполнение. Без модуля значения многочлена степени 300 в случайной точке имеют тысячи цифр, и проверка стоит дороже раскрытия скобок.

Зачем это в работе

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

Сравнение больших данных на расстоянии. Два сервера с копиями файла в гигабайт не пересылают файл: каждый читает его как коэффициенты многочлена, вычисляет значение в общей случайной точке по простому модулю, много большему длины файла, - порядка 2⁶¹ - и пересылает одно число. Это отпечаток из алгоритма Рабина-Карпа, и оценка та же: d / p.

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

Теперь сам

1. Сколько умножений чисел нужно на одну проверку Фрейвалдса для матриц 500 × 500 и сколько - на прямое перемножение?
2. Сколько независимых проверок нужно, чтобы неверное произведение прошло с вероятностью меньше одной миллиардной? 2³⁰ ≈ 1,07 · 10⁹.
3. Два многочлена степени 1000 сравнивают в случайной точке по простому модулю около 10⁹. Оцени вероятность ошибки.
4. В матрице C ошибочны все элементы одного столбца j, остальные верны. С какой вероятностью одна проверка это пропустит?
5. Почему проверку Фрейвалдса нельзя превратить в алгоритм Лас-Вегас, повторяя её до ответа «равно»?