《Magus Night》
求gcd <= q,且lcm >= p的方案数。
考虑容斥:ans = ans{无限制下的序列方案数} - {gcd > q的方案数} - {gcd <= q 且 lcm < p的方案数}
对于无限制下的序列方案数,因为每个数都能选的是[1,m],要求的是序列中数的乘积和。
所以总数量 = (1 + 2 + ... m) * (1 + 2 + ... m) * ... = (1 + 2 + .. m) ^ n
对于gcd > q的方案数,我们枚举gcd容斥dp一下即可。定义dp[i]表示gcd = i时的方案数,就是个经典的dp计数了。
对于gcd <= q 且 lcm < p的方案数。
这里我们完全可以枚举gcd = d,lcm = t来统计方案数:考虑当前gcd = d,lcm = t。
对t质因子分解,p1 ^ k1 * p2 ^ k2 * ... pn ^ kn。对于每个质因子他的幂次只能选[0,k1],并且对于每个上限k,下限1。都要选到。
因为上限不选到k不满足lcm,下限不选到1不满足gcd。
然后关于这个选的方案数也可以容斥计算:定义cal(i,j) - 在[i,j]中随意选的方案数,这个很容易计算。
ans[1,k] = cal(1,k) - cal(2,k) - cal(1,k - 1) + cal(2,k - 1)
f[i] - 下限 = 1,上限 = i时的方案数.
// Author: levil #include<bits/stdc++.h> using namespace std; typedef long long LL; typedef long double ld; typedef pair<int,int> pii; const int N = 2e5 + 5; const int M = 1e6 + 5; const LL Mod = 998244353; #define INF 1e9 #define dbg(ax) cout << "now this num is " << ax << endl; inline long long ADD(long long x,long long y) { if(x + y < 0) return ((x + y) % Mod + Mod) % Mod; return (x + y) % Mod; } inline long long MUL(long long x,long long y) { if(x * y < 0) return ((x * y) % Mod + Mod) % Mod; return x * y % Mod; } inline long long DEC(long long x,long long y) { if(x - y < 0) return (x - y + Mod) % Mod; return (x - y) % Mod; } bool vis[N]; int prime[N],tot = 0; void init() { for(int i = 2;i < N;++i) { if(!vis[i]) prime[++tot] = i; for(int j = 1;j <= tot && prime[j] * i < N;++j) { vis[prime[j] * i] = 1; if(i % prime[j] == 0) break; } } } LL quick_mi(LL a,LL b) { LL re = 1; while(b) { if(b & 1) re = re * a % Mod; a = a * a % Mod; b >>= 1; } return re; } LL dp[N],f[N]; void solve() { init(); LL inv2 = quick_mi(2,Mod - 2); int n,m,p,q;scanf("%d %d %d %d",&n,&m,&p,&q); LL sum = 1LL * (m + 1) * m % Mod * inv2 % Mod; LL ans = quick_mi(sum,n); for(int i = m;i > q;--i) { int up = m / i; LL cnt = 1LL * (up + 1) * up % Mod * inv2 % Mod; dp[i] = quick_mi(i,n) * quick_mi(cnt,n) % Mod; for(int j = i + i;j <= m;j += i) dp[i] = DEC(dp[i],dp[j]); ans = DEC(ans,dp[i]); } for(int i = 1;i <= m;++i) { f[i] = 1; int x = i; for(int j = 1;prime[j] <= x && x != 1;++j) { if(!vis[x]) break; if(x % prime[j] == 0) { LL ma = 1,pre = 0; while(x % prime[j] == 0) { x /= prime[j]; pre = ADD(pre,ma); ma = MUL(ma,prime[j]); } LL p = ADD(pre,ma); LL tmp = DEC(quick_mi(p,n),quick_mi(pre,n)); tmp = DEC(tmp,quick_mi(DEC(p,1),n)); tmp = ADD(tmp,quick_mi(DEC(pre,1),n)); f[i] = MUL(f[i],tmp); } } if(x != 1) { LL pre = 1,ma = x; LL p = ADD(pre,ma); LL tmp = DEC(quick_mi(p,n),quick_mi(pre,n)); tmp = DEC(tmp,quick_mi(DEC(p,1),n)); tmp = ADD(tmp,quick_mi(DEC(pre,1),n)); f[i] = MUL(f[i],tmp); } } for(int i = 1;i <= q;++i) { for(int j = i;j < p;j += i) ans = DEC(ans,quick_mi(i,n) * f[j / i] % Mod); } printf("%lld\n",ans); } int main() { //int ca;scanf("%d",&ca); //while(ca--) { solve(); //} system("pause"); return 0; }