对于一个素数p,先求n=p-1并将其分解为x个素数因子,对于每一个因子q及其指数c,应用Pohlig-Helliman算法求解(a0,a1,a3............ac-1)
根据a=(for i=0 to c-1:ai*q^i)+s*q^(s为某一整数) 求得 x个同余式,最后利用中国剩余定理求解离散对数的解
# -*- coding: utf-8 -*-
_author_ = 'xiao_lu'
import sys
import time
the_a={}
#将n进行分解为素数积
def get_result(n):
i=2
while n!=1:
if n%i==0:
the_a[i]+=1
n=n/i
else: i+=1
#扩展欧几里德算法
def extended_euclid(a,b,x,y):
if b == 0:
x[0] = 1
y[0] = 0
return a
d = extended_euclid(b, a % b, y, x)
y[0] -= a / b * x[0]
return d
#中国剩余定理
def chinese_remainder(b ,w ,len):
ret = 0; n = 1
x = [0]; y = [0]
for i in range(len):
n *= w[i];
for i in range(len):
m = n / w[i];
d = extended_euclid(w[i],m,x,y);
ret = (ret + y[0] * m * b[i]) % n;
return (n + ret % n) % n;
#主函数 Pohilg-helllman算法
if __name__ == '__main__':
t1 = time.time()
print t1
p,a,b = map(int, sys.stdin.readline().split())
n = int(p)-1
#初始化
for i in range(2,10000):
the_a.setdefault(i,0)
# 将n进行分解为素数积
get_result(int(n))
chu=[]; yu=[]
for key, value in the_a.items():
if value != 0:
#Pohlig-Helliman算法
j = 0
A = []; B = []
B.append(int(b))
key1 = key
while j <= int(value) - 1:
x = (int(B[j])** (n/(key1)))%int(p)
for i in range(0, 10000):
if x == (int(a)**((i * n) / int(key)))%int(p):
A.append(i)
break
B.append(( B[j] * (int(a)**((-i) * (int(key)**j))))%int(p))
j += 1
key1 *= key
s=0
for ii in range(0,int(value)):
#求a
s+=int(A[ii])*(int(key)**ii)
chu.append(int(key)**int(value))
yu.append(s)
#利用中国剩余定理求解同余式组
print "log%d %d = %d"%(a,b,(chinese_remainder(yu,chu,len(yu))))
t2 = time.time()
print t2