【GDOI2016模拟3.16】幂（容斥 + 模型复杂转化）

【GDOI2016模拟3.16】幂

(Xin[1,A],Yin[1,B])，问：(x^y)的不用取值个数.
(A,B)都是(10^9)级别.
然后我们开搞.
首先，假设一个合法的(x)可以表示为(x=prod p_i^{q_i})，那么令(d=gcd(q_1,q_2...q_k))
假设(d>1)，显然我们不需要单独考虑，因为它可以继续化简，我们找到最简的那个数然后去一次性处理.
那么此时所有情况都变成了(d=1).
此时再分两种情况讨论，因为我们现在实际上把问题转化为了：

(xin[1,A])，(yin [1,B])，求((c^x)^y)的不同个数.
然后这个问题显然当(x>sqrt{A})时贡献就是(B)。
所以题目转化为：

在(xin[1,sqrt{A}],yin[1,B])时，((c^x)^y)的不同个数.
于是我们可以正难则反，计算一下不合法的个数，假设(k={log_x}^n)，那么问题又可以转化为：

(xin[1,k],yin[1,B])，(xy)的相同个数.
于是这个问题值得思考一下.
咱有一个经典套路：分块计算；
即把总范围([1,KB])变成(K)块，每一块都长的像这样的形式：

[[(i-1)B+1,iB] | (1le i le K) ]

那么现在单独考虑每一块的贡献.
然后假设我们现在处理第(i)个块，那么如果其中某一个数可以被表示为(x*y)的形式.
则必定存在一个(din [i,k])，满足(d|x*y).
因为这样子就必定能满足(xin [1,K],yin [1,B])这两个条件.
所以现在题目再次被转化：

等价于我们要求在区间([(i-1)B+1,iB])范围内有多少个数是([i,k])区间中某一个数的倍数.
这次转化以后我们发现问题就变得熟悉了...
因为(B)比较大，(k)貌似小，所以,,,, 这TM容斥一下不就好了？？
但是，我们发现，(k)最大可能为(30)，这意味着，(2^{30})是接受不了的...
但是我们很容易想到，如果当乘上一个数的lcm没有变化之后，就不需要去弄.
这样子的话时间一下子就降了下来，反正怎么打都能过~

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>

#define ll long long
#define F(i, a, b) for (ll i = a; i <= b; i ++)
#define sqr(x) ((x) * (x))

const ll N = 5e5;

ll ans, sum;

ll A, B, up, P[N];

bool vis[N];

using namespace std;

void Init() {
	#ifndef ONLINE_JUDGE
		freopen("data.in", "r", stdin);
	#endif
	
	scanf("%lld%lld", &A, &B);
	
	up = int(sqrt(A)) + 1;
	
	ans = 1LL * (A - 1) * B + 1;
}

void GetPrime() {
	F(i, 2, N - 1) {
		if (!vis[i])
			P[++ P[0]] = i;
		F(j, 1, P[0]) {
			if (P[j] * i >= N) break;
			vis[P[j] * i] = 1;
			if (i % P[j] == 0) break;
		}
	}
}

ll gcd(ll x, ll y){
	return y == 0 ? x : gcd(y, x % y);
}

void dfs(ll x, ll v, ll st, ll en, ll g, ll k) {
	if (x > k) {
		sum += 1LL * (en / v - st / v) * g;
		return;
	}
	if (x != gcd(v, x)) {
		dfs(x + 1, v * x / gcd(v, x), st, en, g * (- 1), k);
		dfs(x + 1, v, st, en, g, k);
	}
}

ll Solve(ll x) {
	ll cnt = -1;
	for (ll tmp = 1; tmp <= A; tmp *= x, cnt ++);

	sum = 0;
	F(i, 2, cnt)
		dfs(i, 1, (i - 1) * B, i * B, 1, cnt);

	return sum;
}

void Doit() {
	F(i, 2, up) {
		ll tmp = i, yes = 0;
		F(j, 1, P[0]) {
			if (sqr(P[j]) > i) break;
			
			ll cnt = 0;
			while (tmp % P[j] == 0) tmp /= P[j], cnt ++;
			
			yes = gcd(yes, cnt);
			if (yes == 1) break;
		}
		
		if (tmp > 1 || yes == 1)
			ans -= Solve(i);
	}
	
	printf("%lld
", ans);
}

int main() {
	Init();
	
	GetPrime();
	
	Doit();
}