P7486 「StOI2031」彩虹
cyn2006 爆切这题之后疯狂安利我做这题。
stO cyn2006 Orz
这是一种复杂度劣于 std 但是能过的做法。
显然容斥一下转为计算下式
\[\prod_{i=1}^{n}\prod_{j=1}^{m}\operatorname{lcm}(i,j)^{\operatorname{lcm}(i,j)}\\=\prod_{i=1}^{n}\prod_{j=1}^{m}\prod_{d=1}^{n}[\gcd(i,j)=d]\left(\dfrac{ij}{d}\right)^{\frac{ij}{d}}\\=\dfrac{\prod_{i=1}^{n}\prod_{j=1}^{m}\prod_{d=1}^{n}[\gcd(i,j)=d](ij)^{\frac{ij}{d}}}{\prod_{i=1}^{n}\prod_{j=1}^{m}\prod_{d=1}^{n}[\gcd(i,j)=d]d^{\frac{ij}{d}}}
\]
看分母清新亿点就先推分母
【注】定义 \(s1(x)=\sum_{i=1}^{x}i\)，以下不再赘述。
分母
\[\prod_{i=1}^{n}\prod_{j=1}^{m}\prod_{d=1}^{n}d^{[\gcd(i,j)=d]\frac{ij}{d}}\\
=\prod_{d=1}^{n}\prod_{i=1}^{\frac{n}{d}}\prod_{j=1}^{\frac{m}{d}}d^{[\gcd(i,j)=1]ijd}\\
=\prod_{d=1}^{n}d^{d\sum_{i=1}^{\frac{n}{d}}\sum_{j=1}^{\frac{m}{d}}[\gcd(i,j)=1]ij}\\
=\prod_{d=1}^{n}d^{df(\frac{n}{d},\frac{m}{d})}\\
\]
其中
\[f(n,m)=\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=1]ij\\=\sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{x|i,x|j}\mu(x)ij\\=\sum_{x=1}^{n}\mu(x)x^2\sum_{i=1}^{\frac{n}{x}}\sum_{j=1}^{\frac{m}{x}}ij\\=\sum_{x=1}^{n}\mu(x)x^2s1(\frac{n}{x})s1(\frac{m}{x})
\]
对指数整除分块，预处理 \(d^d\) 的前缀积，由于 \(f(n,m)\) 的计算还有一层整除分块，所以这部分是 \(O(n^{\frac{3}{4}})\)。
注意 \(f(n,m)\) 作为指数，是在 \(\bmod mod - 1\) 意义下进行的。
分子
\[\prod_{i=1}^{n}\prod_{j=1}^{m}\prod_{d=1}^{n}(ij)^{[\gcd(i,j)=d]\frac{ij}{d}}\\
=\prod_{i=1}^{n}\prod_{j=1}^{m}\prod_{d=1}^{n}(ijd^2)^{[\gcd(i,j)=1]ijd}\\
=\prod_{d=1}^{n}\prod_{i=1}^{\frac{n}{d}}\prod_{j=1}^{\frac{m}{d}}\left(ij\right)^{[\gcd(i,j)=1]ijd}d^{[\gcd(i,j)=1]2ijd}\\
\]
分两部分算，然后乘起来。
Part1
\[=\prod_{d=1}^{n}\prod_{i=1}^{\frac{n}{d}}\prod_{j=1}^{\frac{m}{d}}[\gcd(i,j)=1]d^{2dij}\\
\]
会发现这个就是分母的平方。
实现的时候可以直接计算分母乘第二部分。
Part2
\[=\prod_{d=1}^{n}\prod_{i=1}^{\frac{n}{d}}\prod_{j=1}^{\frac{m}{d}}\left(ij\right)^{[\gcd(i,j)=1]ijd}\\
\]
指数上那个 \(d\) 暂且不管，反正最后 \(d\) 次方就行。
\[g(n,m)=\prod_{i=1}^{n}\prod_{j=1}^{m}(ij)^{[\gcd(i,j)=1]ij}
\]
再观察一下发现这个东西还是可以拆分的。
\[h(n,m)=\prod_{i=1}^{n}i^{i\sum_{j=1}^{m}[\gcd(i,j)=1]j}
\]
那么
\[g(n,m)=h(n,m)\times h(m,n)
\]
\[\sum_{j=1}^{m}[\gcd(i,j)=1]j=\sum_{j=1}^{m}\sum_{t|i,t|j}\mu(t)j=\sum_{t|i}\mu(t)\sum_{i=1}^{m}[t|i]i=\sum_{t|i}\mu(t)t\sum_{i=1}^{\frac{m}{t}}i=\sum_{t|i}t\mu(t)s1(\frac{m}{t})
\]
\[h(n,m)=\prod_{i=1}^{n}i^{i\sum_{t|i}t\mu(t)s1(\frac{m}{t})}\\=\prod_{t=1}^{n}(\prod_{t|i}i^i)^{s1(\frac{m}{t})t\mu(t)}
\]
预处理 \(R(x,n)=(\prod_{i=1}^{n}(ix)^{ix})^{x\mu(x)}\)，可以发现只有 \(O(n\ln n)\) 种 \((x,n)\)，动态分配内存即可。
\[=\prod_{t=1}^{n}R(t,\dfrac{n}{t})^{s1(\frac{m}{t})}
\]
发现要求出
\[Q(n,m)=\prod_{i=1}^{m}R(i,n)
\]
就可以整除分块了。\(Q\) 仍然是只有 \(O(n\ln n)\) 种。
由于还要快速幂这玩意是 \(O(n^{\frac{3}{4}}\log n)\) 的。
写完发现空间超了 20MB 左右。于是考虑只开 \(Q\) 不开 \(R\)，递推的时候枚举 \(n\)，对所有 \(R(x,n)\) 开个数组维护第一维为 \(x\) 的时候的值，每次 \(n\) 增加的时候都可以递推出 \(R(x,n)\)。
注意到上面递推 \(R\) 的时候有 \(O(n\ln n)\) 次快速幂，我没啥好的处理方法，就加上 \(\mu(x)=0\) 直接设为 \(1\) 的剪枝，写了一发跑起来不是特别慢。
总复杂度是 \(O(tn^{\frac{3}{4}}\log n)\) 还要加个玄学预处理，跑了两秒不到一点，被 std 吊打了，orz 出题人！
#include<bits/stdc++.h>
using namespace std;
#define fi first
#define se second
#define mkp make_pair
#define pb push_back
#define sz(v) (int)(v).size()
typedef long long LL;
typedef double db;
template<class T>bool ckmax(T&x,T y){return x<y?x=y,1:0;}
template<class T>bool ckmin(T&x,T y){return x>y?x=y,1:0;}
#define rep(i,x,y) for(int i=x,i##end=y;i<=i##end;++i)
#define per(i,x,y) for(int i=x,i##end=y;i>=i##end;--i)
inline int read(){
int x=0,f=1;char ch=getchar();
while(!isdigit(ch)){if(ch=='-')f=0;ch=getchar();}
while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
return f?x:-x;
}
const int N = 1000005;
const int mod = 32465177;
inline int qpow(int n, int k) {
int res = 1;
for(; k; k >>= 1, n = (LL)n * n % mod)
if(k & 1) res = (LL)res * n % mod;
return res;
}
/*
f1_k = \prod_{i=1}^{k}i^i
g1_k = (f1_k)^{-1}
f2_k = k^k
g2_k = f2_k^{-1}
f3_k = \sum_{i=1}^{k}\mu(i)i^2
f(n, m) = \sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=1]ij
smu_k = \sum_{i=1}^{k}\mu_i
*/
int T, n;
int f1[N], g1[N], f2[N], g2[N], f3[N], pri[N], pct, mu[N], smu[N];
int *Q[N], pool[N * 20], *mem = pool, R[N];
bool vis[N];
inline void init(const int &n) {
f1[0] = 1, g1[0] = 1;
for(int i = 1; i <= n; ++i)
f2[i] = qpow(i, i), g2[i] = qpow(f2[i], mod - 2),
f1[i] = (LL)f1[i - 1] * f2[i] % mod,
g1[i] = (LL)g1[i - 1] * g2[i] % mod;
mu[1] = 1;
for(int i = 2; i <= n; ++i) {
if(!vis[i]) pri[++pct] = i, mu[i] = -1;
for(int j = 1; j <= pct && i * pri[j] <= n; ++j) {
vis[i * pri[j]] = 1;
if(i % pri[j] == 0) {
mu[i * pri[j]] = 0;
break;
}
mu[i * pri[j]] = -mu[i];
}
}
for(int i = 1; i <= n; ++i)
smu[i] = smu[i - 1] + mu[i],
f3[i] = (f3[i - 1] + (LL)i * i * mu[i] % (mod - 1) + (mod - 1)) % (mod - 1);
for(int i = 1; i <= n; ++i)
Q[i] = mem, mem += n / i + 2, R[i] = 1, Q[i][0] = 1;
for(int i = 1; i <= n; ++i) {
for(int j = 1; j <= n / i; ++j) {
if(mu[j] == 0) R[j] = 1;
else if(mu[j] == 1) R[j] = (LL)R[j] * qpow(f2[i * j], j) % mod;
else R[j] = (LL)R[j] * qpow(g2[i * j], j) % mod;
Q[i][j] = (LL)Q[i][j - 1] * R[j] % mod;
}
}
}
inline int s1(int n, int P = mod) {
return (LL)n * (n + 1) / 2 % P;
}
inline int f(int n, int m) {
int res = 0;
for(int l = 1, r; l <= n; l = r + 1)
r = min(n / (n / l), m / (m / l)),
res = (res + (LL)(f3[r] - f3[l - 1] + mod - 1) *
s1(n / l, mod - 1) % (mod - 1) * s1(m / l, mod - 1) % (mod - 1)) % (mod - 1);
return res;
}
inline int calc1(int n, int m) {
int res = 1;
for(int l = 1, r; l <= n; l = r + 1)
r = min(n / (n / l), m / (m / l)),
res = (LL)res * qpow((LL)f1[r] * g1[l - 1] % mod, f(n / l, m / l)) % mod;
return res;
}
inline int h(int n, int m) {
int res = 1;
for(int l = 1, r; l <= min(n, m); l = r + 1)
r = min(n / (n / l), m / (m / l)),
res = (LL)res *
qpow((LL)Q[n / l][r] * qpow(Q[n / l][l - 1], mod - 2) % mod,
s1(m / l, mod - 1)) % mod;
return res;
}
inline int calc2(int n, int m) {
int res = 1;
for(int l = 1, r; l <= n; l = r + 1)
r = min(n / (n / l), m / (m / l)),
res = (LL)res * qpow((LL)h(n / l, m / l) * h(m / l, n / l) % mod,
(LL)(l + r) * (r - l + 1) / 2 % (mod - 1)) % mod;
return res;
}
inline int solve(int n, int m) {
return (LL)calc1(n, m) * calc2(n, m) % mod;
}
void Main() {
int l = read(), r = read();
printf("%lld\n", (LL)solve(r, r) * solve(l - 1, l - 1) % mod *
qpow(solve(l - 1, r), mod - 3) % mod);
}
signed main() {
for(T = read(), n = read(), init(n); T--; ) Main();
}
路漫漫其修远兮,吾将上下而求索
