Tenka1 Programmer Contest: ModularPowerEquation!! (AC解法(解法なんだからACなのは当たり前だろ))

http://tenka1-2017.contest.atcoder.jp/tasks/tenka1_2017_f

http://omochan.hatenablog.com/entry/2017/10/20/004232の続きです。

普通に蟻本通りにやればいいだけでした…。
x≡k(modm)
x≡k'(modm')を満たすxは
x=mt+kとして
mt≡k'-k(modm')をtについて解けばいいだけでした。k,m,k',m'がすべて10^9より小さければ、xはオーバーフローせずに求められます。
具体的にはmt'≡g(modm')(g=gcd(m,m'))となるt'(modm')を求めてからt=(k'-k)/g*t'とすればいいだけです。

modの基本定理のおさらい。

  • g|Mのときga≡gb(modM)<=>a≡b(modM/g)
  • gとMが互いに素の時,ga≡gb(modM)<=>a≡b(modM)

二つとも両辺に寄せてからg(a-b)=Mtとかすれば導出できます。これが一番modの中で基礎となる式かもしれない。modの積と商の変形ルールを述べています。
ちなみに和(と差)は当たり前ですが、
a≡b(modM)<=>a+c≡b+c(modM)です。

四則演算の変形+中国剰余定理でだいぶmodは使いやすくなった?

ll phi(ll M) {
	ll A = M;
	for(ll i = 2; i * i <= A; i++) {
		if(A % i == 0) {
			M = M / i * (i - 1);
			while(A % i == 0) A /= i;
		}
	}
	if(A != 1) M = M / A * (A - 1);
	return M;
}

ll gcd(ll a, ll b) {
	return b == 0 ? a : gcd(b, a % b);
}

ll extgcd(ll a, ll b, ll& x, ll& y) { //get(ax + by = gcd(a, b))
	if(b == 0) {
		x = 1; y = 0;
		return a;
	}
	ll res = extgcd(b, a % b, y, x);
	y -= (a / b) * x;
	return res;
}

ll mod_pow(ll a, ll n, ll mod) { //get a^n mod
	if(n == 0) return 1;
	ll res = mod_pow(a, n / 2, mod);
	if(n % 2 == 0) return res * res % mod;
	else return a * res % mod * res % mod;
}

ll mod_inv(ll a, ll m) { //get x such that ax=gcd(a, m)(modm)
	ll x, y;
	extgcd(a, m, x, y);
	return (x % m + m) % m;
}

int Q;
ll A, M;

ll loop(ll A, ll M) {
	if(M == 1) return 10000;
	ll K = loop(A, phi(M));
	ll tK = mod_pow(A, K, M);
	return tK + (tK < 100 ? M : 0);
}

void solve() {
	cin >> Q;
	while(Q--) {
		cin >> A >> M;
		ll pM = phi(M);
		ll K = loop(A, M) % M;
		ll pK = loop(A, pM) % pM;
		ll pt = mod_inv(M, pM);
		ll g = gcd(M, pM), l = M * pM / g;
		ll t = ((pK - K) / g * pt) % pM;
		//debug(M, pM, K, pK, pt, g, l, t);
		ll pans = ((t * M % l + K) % l + l) % l;
		ll ans = pans < 100 ? pans + l : pans;
		cout << ans << "\n";
		//debug(mod_pow(A, ans, M), ans % M);
	}
}