Discrete Logarithm

Definition

Prerequisite knowledge: order and primitive roots.

Discrete logarithms are defined in a similar way to logarithms. Take the positive integer modulus $m$ with primitive roots, and let one of its primitive roots be $g$. For the integer $a$ that satisfies $(a,m)=1$, we know that there must be a unique integer $0\leq k<\varphi(m)$ such that

We call this $k$ the discrete logarithm of $g$ modulo $m$, denoted as $k=\operatorname{ind}_g a$, or $\operatorname{ind} a$ without causing confusion.

Apparently $\operatorname{ind}_g 1=0$, $\operatorname{ind}_g g=1$.

Properties

The properties of discrete logarithms also have many similarities with logarithms.

Baby-Step Giant-Step Algorithm

At present, there is still no polynomial time classical algorithm for the discrete logarithm problem (the input size of the discrete logarithm problem is the number of digits of the input data). In cryptography, based on this, people have designed many asymmetric encryption algorithms, such as Ed25519.

In algorithm competitions, BSGS (baby-step giant-step algorithm) is often used to solve discrete logarithm problems. Formally speaking, for $a,b,m\in\mathbf{Z}^+$, the algorithm can be solved in $O(\sqrt{m})$ time

Here $a\perp m$. The solution $x$ satisfies $0 \le x < m$. (Note that $m$ is not necessarily a prime number)

Algorithm Description

Let $x = A \left \lceil \sqrt m \right \rceil - B$, among them $0\le A,B \le \left \lceil \sqrt m \right \rceil$, then we have $a^{A\left \lceil \sqrt m \right \rceil -B} \equiv b \pmod m$, and with a slight transformation, we have $a^{A\left \lceil \sqrt m \right \rceil} \equiv ba^B \pmod m$.

What we know is $a,b$, so we can first calculate all the values ​​of $ba^B$ on the right side of the equation, enumerate $B$, save it with hash/map, then calculate $a^{A\left \lceil \sqrt m \right \rceil}$ one by one, enumerate $A$, and find whether there is an equal value $ba^B$, so we can get all $x$, $x=A \left \lceil \sqrt m \right \rceil - B$.

Note that $A,B$ is smaller than $\left \lceil \sqrt m \right \rceil$, so the time complexity is $\Theta\left (\sqrt m\right )$, and using map will result in one more $\log$.

Why does > [!note]- require that $a$ and $m$ are relatively prime?

Note that what we find is $A,B$. We need to ensure that $a^{A\left \lceil \sqrt m \right \rceil} \equiv ba^B \pmod m$ can be pushed back to $a^{A\left \lceil \sqrt m \right \rceil -B} \equiv b \pmod m$. The latter formula is obtained by dividing the left and right sides of the front formula by $a^B$, so there must be $a^B \perp m$, that is, $a\perp m$.

Extended BSGS Algorithm

For $a,b,m\in\mathbf{Z}^+$, solve

Here $a,m$ are not necessarily coprime.

When $(a, m)=1$, $a$ has an inverse element modulo $m$, so it can be solved using the BSGS algorithm. So we tried to find ways to make them mutually prime.

Specifically, let $d_1=(a, m)$. If $d_1\nmid b$, then the original equation has no solution. Otherwise we divide the equation by $d_1$ at the same time and get

If $a$ and $\frac{m}{d_1}$ are still not mutually prime, divide again, let $d_2=\left(a, \frac{m}{d_1}\right)$. If $d_2\nmid \frac{b}{d_1}$, the equation has no solution; otherwise divide by $d_2$ at the same time to get

In the same way, continue to judge in this way until $a\perp \dfrac{m}{d_1d_2\cdots d_k}$.

Remember $D=\prod_{i=1}^kd_i$, so the equation becomes like this:

Since $a\perp\dfrac{m}{D}$, $\dfrac{a^k}{D}\perp \dfrac{m}{D}$ follows. In this way, $\dfrac{a^k}{D}$ has an inverse element, so move it to the right-hand side. This is an ordinary BSGS problem. So after solving $x-k$ and adding $k$, it is the solution of the original equation.

Note that it is not excluded that the solution is less than or equal to $k$, so before eliminating factors, do the $\Theta(k)$ enumeration and directly verify $a^i\equiv b \pmod m$, so that this situation can be avoided.

Fast Discrete Logarithm with Range Preprocessing

The time complexity of the previous BSGS algorithm is a single $O(\sqrt m)$, and the efficiency is low when the query magnitude is large. If the modulus solved each time is a fixed prime number $p$, we have a fast algorithm based on range preprocessing.

We already know $\operatorname{ind}_g(ab)\equiv\operatorname{ind}_g a+\operatorname{ind}_g b\pmod{p-1}$, so we can only calculate the discrete logarithm of all prime numbers through the BSGS algorithm. The discrete logarithm of the composite number can be converted into the sum of the discrete logarithms of several known prime numbers through this formula. At this time, the complexity is still not optimal. We consider preprocessing only a part of the discrete logarithms. Specifically, we preprocess the discrete logarithms from $1$ to $L = \lfloor\sqrt p\rfloor + 1$. Note that the BSGS block length $B$ at this time cannot be taken $O(\sqrt{L})$, because the complexity of the BSGS preprocessing (insertion into the hash table) part is $O(B)$, and the query requires a total of $O(\pi(L))$ times, so the total time complexity is $O\left(B+\dfrac{\pi(L)p}{B}\right)$, and $B=O(\sqrt{\pi(L)p})$ is optimal at this time. Using the prime number theorem to $\pi(n)\sim\dfrac{n}{\log n}$, the total preprocessing time complexity can be balanced to $O\left(\dfrac{p^{3/4}}{\log^{1/2} p}\right)$.

Next is how to find the answer. Assume that the current requirement is $\operatorname{ind}_g y$, if $y\le L$, return directly, otherwise set $p=vy+r$, then $v=\left\lfloor\dfrac{p}{y}\right\rfloor<L$, $r=p\bmod y$, $y=\dfrac{p-r}{v}$, thus

Note that $\operatorname{ind}_g (p-1)=(p-1)/2$, so we only need to recursively calculate the discrete logarithm of $r$.

We can also consider another expression of $y$, noting $p=vy+r=(v+1)y+r-y$, then $y=\dfrac{p-r+y}{v+1}$, thus

We have $v+1 \le L$, so we only need to recursively calculate the discrete logarithm of $y-r$.

Combining these two calculation methods, we have $\min\{r,y-r\}\le \dfrac{y}{2}$, so the smaller recursive calculation can achieve the query complexity of $O(\log p)$.

So far we have obtained an algorithm with a time complexity of $O\left(\dfrac{p^{3/4}}{\log^{1/2} p}\right)-O(\log p)$.

#include <cmath>
#include <iostream>
#include <unordered_map>
using namespace std;

const int N = 1'000'000;  // sqrt(P * sqrt(P) / ln(P))

unordered_map<int, int> M;
int Lg[N], p[N], t;
bool vis[N];
int P, g, B, sq, LP_1, inv;

int fadd(int x, int y, int P) {
  x += y;
  if (x >= P) x -= P;
  return x;
}

int fsub(int x, int y, int P) {
  x -= y;
  if (x < 0) x += P;
  return x;
}

int fmul(int x, int y, int P) { return 1ll * x * y % P; }

int qpow(int x, int y) {
  int ans = 1;
  while (y) {
    if (y & 1) ans = fmul(ans, x, P);
    x = fmul(x, x, P);
    y >>= 1;
  }
  return ans;
}

int calc(int x) {
  int s = x;
  for (int i = 0; i <= P / B; ++i) {
    if (M.find(s) != M.end()) return i * B + M[s];
    s = fmul(s, inv, P);
  }
  return -1;
}

void init() {
  int s = 1;
  for (int i = 0; i < B; ++i) {
    if (M.find(s) != M.end()) break;
    M[s] = i, s = fmul(s, g, P);
  }
  inv = qpow(qpow(g, B), P - 2);
  for (int i = 2; i <= sq; ++i) {
    if (!vis[i]) {
      p[++t] = i;
      Lg[i] = calc(i);
    }
    for (int j = 1; j <= t && p[j] * i <= sq; ++j) {
      vis[p[j] * i] = true;
      Lg[p[j] * i] = fadd(Lg[p[j]], Lg[i], P - 1);
      if (i % p[j] == 0) break;
    }
  }
}

int solve(int y) {
  if (y <= sq) return Lg[y];
  int v = P / y, r = P % y;
  if (r < y - r) return fadd(fsub(solve(r), Lg[v], P - 1), LP_1, P - 1);
  return fsub(solve(y - r), Lg[v + 1], P - 1);
}

int main() {
  ios::sync_with_stdio(false);
  cin.tie(nullptr);
  cin >> P >> g;
  sq = sqrt(P) + 1;
  B = sqrt(1ll * P * sqrt(P) / log(P));
  init();
  LP_1 = (P - 1) / 2;  // g ^ LP_1 = P - 1 (mod P)
  int T;
  cin >> T;
  while (T--) {
    int x;
    cin >> x;
    cout << solve(x) << '\n';
  }
  return 0;
}

Exercises

• SPOJ MOD template
• SDOI2013 random number generator
• SGU261 Discrete Roots template
• SDOI2011 calculator template
• Luogu4195【template】exBSGS/Spoj3105 Mod template
• Codeforces - Lunar New Year and a Recursive Sequence
• LOJ6542 discrete logarithm index calculus method, non-template

**Part of the content and code of this page are adapted from the blog post Дискретное извлечение корня and its English-language edition Discrete Root. The copyright agreement for the Russian version is Public Domain + Leave a Link; the copyright agreement for the English version is CC-BY-SA 4.0. **

References

1. Discrete logarithm - Wikipedia
2. Pan Chengdong, Pan Chengbiao. Elementary number theory.
3. Feng Keqin． Elementary number theory and its applications.