转自：http://blog.miskcoo.com/2015/05/discrete-logarithm-problem

扩展大步小步法解决离散对数问题

May 16, 2015 miskcoo Algorithm, Math

离散对数（Discrete Logarithm）问题是这样一个问题，它是要求解模方程

a x \equiv b (mod m)

这个问题是否存在多项式算法目前还是未知的，这篇文章先从 m 是质数开始介绍大步小步法（Baby Step Giant Step）来解决它，之后再将其应用到 m 是任意数的情况。这个算法可以在 O(m−−√) 的时间内计算出最小的 x ，或者说明不存在这样一个 x

题目链接：BZOJ-2480、SPOJ-MOD、BZOJ-3239

首先解决 m=p 是质数的情况，我们可以设 x=A⌈p–√⌉+B ，其中 0≤B<⌈p–√⌉ , 0≤A<⌈p–√⌉ ，这样的话就变成求解 A 和 B 了，方程也变成

a A ⌈ p \sqrt ⌉ + B \equiv b (mod p)

可以在两边同时乘以 aB 的逆元，由于 p 是质数，这个逆元一定存在，于是方程变成

a A ⌈ p \sqrt ⌉ \equiv b \cdot a - B (mod p)

由于 A,B 都是 O(p–√) 级别的数，可以先计算出右边这部分的值，存入 Hash 表，然后计算左边的值，在 Hash 表中查找，只要按照从小到大的顺序如果有解就能够找到最小的解，由于两边都只有 O(p–√) 个数，因此时间复杂度是 O(p–√) 的，这样 m 是质数的情况就解决了

一个优化：我们可以设 x=A⌈p–√⌉−B ，其中 0≤B<⌈p–√⌉ , 0<A≤⌈p–√⌉+1 ，这样的话化简后的方程就是

a A ⌈ p \sqrt ⌉ \equiv b \cdot a B (mod p)

就可以不用求出逆元，要注意只是不用求出逆元，而不是没有用到逆元的存在

现在来看 m 不是质数的情况，同样可以设 x=A⌈m−−√⌉+B ，根据上面的推导，会发现需要用到的性质就是 aB 的逆元存在，所以当 m 和 a 互质的时候这个方法仍然有效！

如果 (m,a)≠1 该怎么办呢？我们要想办法把方程转化为 (m,a)=1 的情况

把要求的模方程写成另外一种形式

a x + k m = b, k \in Z

设 g=(a,m) ，这样的话可以确定如果 g∤b 那么该方程一定无解，所以当 g∣b 的时候，在方程左右两边同时除以 g

a g a x - 1 + k m g = b g, k \in Z

这样便消去了一个因子，得到方程

a g a x - 1 \equiv b g (mod m g)

令 m′=mg,b′=bg(ag)−1 （这里不可以把 g 消掉），就可以得到新的方程

a x' \equiv b' (mod m')

得到解之后原方程的解 x=x′+1 ，不断重复这个过程最后一定会得到一个可以解的方程，套用刚刚的大步小步法解出后即可。要注意的是在这个过程中如果某一步发现 b′=1 ，那么就可以直接退出，因为这时候已经得到了解

NOTE：上面这个过程是可能执行多次的比如说 (a,m)=(6,16)→(6,8)→(6,4)→(6,2)

下面的是代码，题目是文章开头给的链接

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
/* BZOJ-2480: Spoj3105 Mod
*   扩展大步小步  */
#include <cstdio>
#include <cmath>
#include <map>
 
int  mod_pow ( long  long  x ,  long  long  p ,  long  long  mod_v )
{
 long  long  v  =  1 ;
 while ( p )
 {
 if ( p  &  1 )  v  =  x  *  v  %  mod_v ;
 x  =  x  *  x  %  mod_v ;
 p  >> =  1 ;
 }
 
 return  v ;
}
 
int  gcd ( int  a ,  int  b )
{
 return  b  ?  gcd ( b ,  a  %  b )  :  a ;
}
 
int  baby_step_giant_step ( int  a ,  int  b ,  int  p )
{
 a  %=  p ,  b  %=  p ;
 if ( b  ==  1 )  return  0 ;
 int  cnt  =  0 ;
 long  long  t  =  1 ;
 for ( int  g  =  gcd ( a ,  p ) ;  g  !=  1 ;  g  =  gcd ( a ,  p ) )
 {
 if ( b  %  g )  return  - 1 ;
 p  /=  g ,  b  /=  g ,  t  =  t  *  a  /  g  %  p ;
 ++ cnt ;
 if ( b  ==  t )  return  cnt ;
 }
 
 std :: map < int ,  int >  hash ;
 int  m  =  int ( sqrt ( 1.0  *  p )  +  1 ) ;
 long  long  base  =  b ;
 for ( int  i  =  0 ;  i  !=  m ;  ++ i )
 {
 hash [ base ]  =  i ;
 base  =  base  *  a  %  p ;
 }
 
 base  =  mod_pow ( a ,  m ,  p ) ;
 long  long  now  =  t ;
 for ( int  i  =  1 ;  i  <=  m  +  1 ;  ++ i )
 {
 now  =  now  *  base  %  p ;
 if ( hash . count ( now ) )
 return  i  *  m  -  hash [ now ]  +  cnt ;
 }
 
 return  - 1 ;
}
 
int  main ( )
{
 int  a ,  b ,  p ;
 while ( std :: scanf ( "%d %d %d" ,  &a ,  &p ,  &b ) ,  p )
 {
 int  ans  =  baby_step_giant_step ( a ,  b ,  p ) ;
 if ( ans  ==  - 1 )  std :: puts ( "No Solution" ) ;
 else  std :: printf ( "%d\n" ,  ans ) ;
 }
 return  0 ;
}