[ZOJ3254] MON 9.2009Secret Code

A^x = D (mod P)

0 <= x <= M, here M is a given integer. 

1 <= A, P < 2^31, 0 <= D < P, 1 <= M < 2^63

----------------------------------------------------------

裸拓展baby step giant step

先转成非拓展版，然后如果在转的过程中就出解了，则return 1,否则就找出=D的起点和循环节长度，然后直接求解。

#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
const int Hash_MOD = 1000003;
typedef long long LL;
typedef unsigned long long ULL;
LL times;
struct link
{
	LL link;
	int val,next;
}es[100000];
struct triple
{
	LL x,y,g;
	triple(const LL _x = 0,const LL _y = 0,const LL _g = 0):
		x(_x),y(_y),g(_g){}
};
int H[Hash_MOD * 2],ec;
triple exgcd(const LL a,const LL b)
{
	if (!b) return triple(1,0,a);
	const triple last(exgcd(b,a%b));
	return triple(last.y,last.x - a / b * last.y,last.g);
}
LL rem_equ(const LL a,const LL b,const LL c)
{
	//ax == c (mod b)
	//ax mod b == c
	const triple tmp(exgcd(a,b));
	const LL MOD = (b / tmp.g);
	return ((tmp.x * (c / tmp.g) % MOD) + MOD) % MOD;
}
LL find(const LL x)
{
	for (int i = H[x%Hash_MOD]; i; i = es[i].next)
		if (es[i].link == x) return es[i].val;
	return -1;
}
void insert(LL x,const LL v)
{
	if (find(x) != -1) return;
	const LL key = x % Hash_MOD;
	es[++ec].next = H[key];
	H[key] = ec;
	es[ec].link = x;
	es[ec].val = v;
}
LL BSGS(LL A,LL P,LL D)
{
	LL AA = 1 % P,x = 1 % P,MOD = P,DD = D,m = static_cast<LL>(std::ceil(std::sqrt((double)P)));
	times = 0;
	for (triple tmp; (tmp = exgcd(A,P)).g != 1; ) {
		if (x == DD) return times++;
		if (D % tmp.g) return -1;
		P /= tmp.g;
		D /= tmp.g;
		(AA *= A / tmp.g) %= P;
		++times;
		(x *= A) %= MOD;
	}
	A %= P;
	ec = 0;
	memset(H,0,sizeof H);
	LL tmp = 1 % P;
	for (int i = 0; i < m; ++i,(tmp *= A) %= P)	
		insert(tmp,i);
	x = 1 % P;
	for (LL i = 0; i < m; ++i,(x *= tmp) %= P) {
		const int j = find(rem_equ(AA*x%P,P,D));
		if (j != -1) return i * m + j + times;
	}
	return -1;
}
LL getphi(LL x)
{
	LL res = 1;
	for (LL i = 2; i * i <= x; ++i)
		if (x % i == 0) {
			x /= i;
			res *= i - 1;
			while (x % i == 0) x /= i,res *= i;
		}
	if (x > 1) res *= x - 1;
	return res;
}
LL power(LL a,LL k,LL MOD)
{
	//a ^ k % MOD
	LL res = 1 % MOD;
	for (; k; k/=2,(a *= a) %= MOD)
		if (k & 1) (res *= a) %= MOD;
	return res;
}
LL calc_len(const LL start,const LL A,const LL P,const LL D)
{
	//A ^ (start + *this) == A ^ start == D (mod P)
	LL phi = getphi(P),res = phi;
	for (LL i = 2; i * i <= phi; ++i) {
		for (; phi % i == 0; phi /= i);
		for (; res % i == 0 && power(A,start + res/i,P) == D; res /= i);
	}
	if (phi > 1)
		for (; res % phi == 0 && power(A,start + res/phi,P) == D; res /= phi);
	return res;
}
ULL work(const LL A,const LL P,const LL D,const ULL M)
{
	LL start = BSGS(A,P,D);
	//printf("%lld\n",start);
	if (start == -1 || start > M) return 0;
	else if (start < times) return 1;
	// LL phi=getphi(P);
	// if (power(A,start + phi,P) != D) return 1;
	ULL len = calc_len(start,A,P,D);
	return (M - start) / len + 1;
}
int main()
{
#ifndef ONLINE_JUDGE
	freopen("3254.in","r",stdin);
	freopen("3254.out","w",stdout);
#endif
	LL A,P,D;
	ULL M;
	while (~scanf("%lld%lld%lld%llu",&A,&P,&D,&M))
		printf("%llu\n",work(A%P,P,D,M));
}