E010. Дискретное логарифмирование

Источник: e-maxx.ru/algo, страница PDF 31.

Problème дискретного логарифмирования заключается в том, чтобы по данным целым

,

,

решить уравнение:

где

и

— взаимно просты (примечание: если они не взаимно просты, то описанный ниже Algorithme является некорректным; хотя, предположительно, его можно модифицировать, чтобы он по-прежнему работал). Здесь описан Algorithme, известный как "baby-step-giant-step algorithm", предложенный

Шэнксом (Shanks) в 1971 г., работающий за время за

. Часто этот Algorithme просто

называют Algorithmeом "meet-in-the-middle" (потому что это одно из классических применений техники "meet-in- the-middle": "разделение задачи пополам").

Algorithme

Итак, мы имеем уравнение:

где

и

взаимно просты.

Преобразуем уравнение. Положим

где

— это заранее выбранная константа (как её выбирать в зависимости от

, мы поймём чуть позже). Иногда

называют "giant step" (поскольку увеличение его на единицу увеличивает

сразу на

), а в противоположность ему

— "baby step".

Очевидно, что любое

(из промежутка

— понятно, что такого диапазона значений будет достаточно)

можно представить в такой форме, причём для этого будет достаточно значений: Тогда уравнение принимает вид:

откуда, пользуясь тем, что

и

взаимно просты, получаем: Чтобы решить исходное уравнение, нужно find соответствующие значения

и

, чтобы значения левой и правой

частей совпали. Иначе говоря, надо решить уравнение: Эта Problème решается с помощью метода meet-in-the-middle следующим образом. Первая фаза Algorithmeа:

посчитаем значения функции

для всех значений аргумента

, и отсортируем эти значения. Вторая фаза

Algorithmeа: будем перебирать значение второй переменной

, вычислять вторую функцию

, и искать это значение

среди предвычисленных значений первой функции с помощью бинарного поиска.

Asymptotic complexity

Сначала оценим время вычисления каждой из функций

и

. И та, и другая содержит возведение в

степень, которое можно выполнять с помощью Algorithmeа бинарного возведения в степень. Тогда обе этих функции

мы можем вычислять за время

.

Сам Algorithme в первой фазе содержит вычисление функции

для каждого возможного значения

и

дальнейшую сортировку значений, что даёт нам асимптотику:

Во второй фазе Algorithmeа происходит вычисление функции

для каждого возможного значения

и бинарный

поиск по tableauу значений

, что даёт нам асимптотику:

Теперь, когда мы сложим эти две асимптотики, у нас получится

, умноженный на сумму

и

, и

практически очевидно, что минимум достигается, когда

, т.е. для оптимальной работы Algorithmeа константу

следует выбирать так: Тогда Asymptotic complexity Algorithmeа принимает вид:

Примечание. Мы могли бы обменять ролями

и

(т.е. на первой фазе вычислять значения функции

, а а второй

—

), однако легко понять, что результат от этого не изменится, и асимптотику этим мы никак не улучшим.

Implémentation

Простейшая Implémentation

Функция

выполняет бинарное возведение числа

в степень

по модулю

, см. Бинарное возведение

в степень.

Функция

производит собственно Solution задачи. Эта функция returns ответ (number в промежутке

), точнее говоря, один из ответов. Функция вернёт

, если решения не существует.

int powmod (int a, int b, int m) {

int res = 1;

while (b > 0)

if (b & 1) {

res = (res * a) % m;

--b;

}

else {

a = (a * a) % m;

b >>= 1;

}

return res % m;

}

int solve (int a, int b, int m) {

int n = (int) sqrt (m + .0) + 1;

map<int,int> vals;

for (int i=n; i>=1; --i)

vals[ powmod (a, i * n, m) ] = i;

for (int i=0; i<=n; ++i) {

int cur = (powmod (a, i, m) * b) % m;

if (vals.count(cur)) {

int ans = vals[cur] * n - i;

if (ans < m)

return ans;

} }

return -1;

} Здесь мы для удобства при реализации первой фазы Algorithmeа воспользовались структурой данных "map" (красно-чёрным arbreм), которая для каждого значения функции

хранит аргумент

, при котором

это значение достигалось. При этом если одно и то же значение достигалось несколько раз, записывается наименьший из всех аргументов. Это сделано для того, чтобы впоследствии, на второй фазе Algorithmeа, нашёлся ответ в

промежутке

.

given, что аргумент функции

на первой фазе у нас перебирался от единицы и до

, а аргумент функции

на второй фазе перебирается от нуля до

, то в итоге мы покрываем всё множество возможных ответов, т.к.

отрезок

содержит в себе промежуток

. При этом отрицательным ответ получиться не мог, а

ответы, большие либо равные

мы можем игнорировать — всё равно должны находиться соответствующие им

ответы из промежутка

. Эту функцию можно изменить на тот случай, если it is required находить все решения задачи дискретного логарифма. Для этого надо заменить "map" на какую-либо другую структуру данных, позволяющую хранить для одного аргумента сразу несколько значений (наExemple, "multimap"), и соответствующим образом изменить код второй фазы.

Улучшенная Implémentation

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

деле, достаточно один раз посчитать величину

, и потом просто домножать на неё. Таким образом, логарифм в асимптотике по-прежнему останется, но это будет только логарифм, связанный со

структурой данных

(т.е., в терминах Algorithmeа, с сортировкой и бинарным поиском значений) — т.е. это

будет логарифм от

, что на практике даёт заметное ускорение.

int solve (int a, int b, int m) {

int n = (int) sqrt (m + .0) + 1;

int an = 1;

for (int i=0; i<n; ++i)

an = (an * a) % m;

map<int,int> vals;

for (int i=1, cur=an; i<=n; ++i) {

if (!vals.count(cur))

vals[cur] = i;

cur = (cur * an) % m;

}

for (int i=0, cur=b; i<=n; ++i) {

if (vals.count(cur)) {

int ans = vals[cur] * n - i;

if (ans < m)

return ans;

}

cur = (cur * a) % m;

}

return -1;

}

Наконец, если модуль

достаточно мал, то можно и вовсе избавиться от логарифма в асимптотике — просто

заведя вместо

обычный tableau. Также можно вспомнить про хеш-таблицы: в среднем они работают также за

, что в целом даёт

асимптотику

.

using System;
using System.Collections.Generic;
using System.Linq;

public static class AlgorithmDraft
{
    // Auto-generated C# draft from the original e-maxx C/C++ listing. Review before production use.
    int powmod (int a, int b, int m) {
            int res = 1;
            while (b > 0)
                    if (b & 1) {
                            res = (res * a) % m;
                            --b;
                    }
                    else {
                            a = (a * a) % m;
                            b >>= 1;
                    }
            return res % m;
    }
    int solve (int a, int b, int m) {
            int n = (int) sqrt (m + .0) + 1;
            map<int,int> vals;
            for (int i=n; i>=1; --i)
                    vals[ powmod (a, i * n, m) ] = i;
            for (int i=0; i<=n; ++i) {
                    int cur = (powmod (a, i, m) * b) % m;
                    if (vals.count(cur)) {
                            int ans = vals[cur] * n - i;
                            if (ans < m)
                                    return ans;
                    }
            }
            return -1;
    }
    int solve (int a, int b, int m) {
            int n = (int) sqrt (m + .0) + 1;
            int an = 1;
            for (int i=0; i<n; ++i)
                    an = (an * a) % m;
            map<int,int> vals;
            for (int i=1, cur=an; i<=n; ++i) {
                    if (!vals.count(cur))
                            vals[cur] = i;
                    cur = (cur * an) % m;
            }
            for (int i=0, cur=b; i<=n; ++i) {
                    if (vals.count(cur)) {
                            int ans = vals[cur] * n - i;
                            if (ans < m)
                                    return ans;
                    }
                    cur = (cur * a) % m;
            }
            return -1;
    }
}

int powmod (int a, int b, int m) {
        int res = 1;
        while (b > 0)
                if (b & 1) {
                        res = (res * a) % m;
                        --b;
                }
                else {
                        a = (a * a) % m;
                        b >>= 1;
                }
        return res % m;
}
int solve (int a, int b, int m) {
        int n = (int) sqrt (m + .0) + 1;
        map<int,int> vals;
        for (int i=n; i>=1; --i)
                vals[ powmod (a, i * n, m) ] = i;
        for (int i=0; i<=n; ++i) {
                int cur = (powmod (a, i, m) * b) % m;
                if (vals.count(cur)) {
                        int ans = vals[cur] * n - i;
                        if (ans < m)
                                return ans;
                }
        }
        return -1;
}
int solve (int a, int b, int m) {
        int n = (int) sqrt (m + .0) + 1;
        int an = 1;
        for (int i=0; i<n; ++i)
                an = (an * a) % m;
        map<int,int> vals;
        for (int i=1, cur=an; i<=n; ++i) {
                if (!vals.count(cur))
                        vals[cur] = i;
                cur = (cur * an) % m;
        }
        for (int i=0, cur=b; i<=n; ++i) {
                if (vals.count(cur)) {
                        int ans = vals[cur] * n - i;
                        if (ans < m)
                                return ans;
                }
                cur = (cur * a) % m;
        }
        return -1;
}

import java.util.*;
import java.math.*;

public class AlgorithmDraft {
    // Auto-generated Java draft from the original e-maxx C/C++ listing. Review before production use.
    int powmod (int a, int b, int m) {
            int res = 1;
            while (b > 0)
                    if (b & 1) {
                            res = (res * a) % m;
                            --b;
                    }
                    else {
                            a = (a * a) % m;
                            b >>= 1;
                    }
            return res % m;
    }
    int solve (int a, int b, int m) {
            int n = (int) sqrt (m + .0) + 1;
            map<int,int> vals;
            for (int i=n; i>=1; --i)
                    vals[ powmod (a, i * n, m) ] = i;
            for (int i=0; i<=n; ++i) {
                    int cur = (powmod (a, i, m) * b) % m;
                    if (vals.count(cur)) {
                            int ans = vals[cur] * n - i;
                            if (ans < m)
                                    return ans;
                    }
            }
            return -1;
    }
    int solve (int a, int b, int m) {
            int n = (int) sqrt (m + .0) + 1;
            int an = 1;
            for (int i=0; i<n; ++i)
                    an = (an * a) % m;
            map<int,int> vals;
            for (int i=1, cur=an; i<=n; ++i) {
                    if (!vals.count(cur))
                            vals[cur] = i;
                    cur = (cur * an) % m;
            }
            for (int i=0, cur=b; i<=n; ++i) {
                    if (vals.count(cur)) {
                            int ans = vals[cur] * n - i;
                            if (ans < m)
                                    return ans;
                    }
                    cur = (cur * a) % m;
            }
            return -1;
    }
}