【BZOJ】3529: [Sdoi2014]数表
题意：求
$$\sum_{i=1}^{n} \sum_{j=1}^{m} \sum_{d|(i, j)} d 且 (\sum_{d|(i, j)} d)<=a$$
n, m<=1e5，q次询问，q<=2*1e4
#include <bits/stdc++.h>
using namespace std;
const int N=1e5+10, MN=1e5, YU=(1u<<31)-1;
int c[N], mx;
void upd(int x, int s) { for(; x<=mx; x+=x&-x) c[x]+=s; }
int sum(int x) { int r=0; for(; x; x-=x&-x) r+=c[x]; return r; }
bool np[N];
int p[N], pcnt, mu[N], f[N], last[N];
void init() {
mu[1]=f[1]=1;
int i, j;
for(i=2; i<=MN; ++i) {
if(!np[i]) p[++pcnt]=i, mu[i]=-1, f[i]=1+i, last[i]=1;
for(j=1; j<=pcnt; ++j) {
int t=p[j]*i; if(t>MN) break;
np[t]=1;
if(i%p[j]==0) { mu[t]=0; f[t]=last[i]+p[j]*f[i]; last[t]=last[i]; break; }
mu[t]=-mu[i]; f[t]=f[i]*(1+p[j]);
last[t]=f[i];
}
}
}
void update(int d, int s) {
for(int i=d, j=1; i<=mx; i+=d, ++j) upd(i, s*mu[j]);
}
struct A { int n, m, a, id; }q[N];
struct B { int f, id; }F[N];
inline bool cmp1(const A &a, const A &b) { return a.a<b.a; }
inline bool cmp2(const B &a, const B &b) { return a.f<b.f; }
int Ans[N];
int main() {
init();
int Q; scanf("%d", &Q);
for(int i=1; i<=Q; ++i) scanf("%d %d %d", &q[i].n, &q[i].m, &q[i].a), q[i].id=i, mx=max(mx, max(q[i].n, q[i].m));
sort(q+1, q+1+Q, cmp1);
for(int i=1; i<=mx; ++i) F[i].id=i, F[i].f=f[i];
sort(F+1, F+1+mx, cmp2);
int n, m, a, now=1, pos, ans;
for(int k=1; k<=Q; ++k) {
n=q[k].n, m=q[k].m, a=q[k].a; if(n>m) swap(n, m);
while(now<=mx && F[now].f<=a) update(F[now].id, F[now].f), ++now;
pos=ans=0;
for(int i=1; i<=n; i=pos+1) {
pos=min(n/(n/i), m/(m/i));
ans+=(sum(pos)-sum(i-1))*(n/i)*(m/i);
}
Ans[q[k].id]=ans&YU;
}
for(int i=1; i<=Q; ++i) printf("%d\n", Ans[i]);
return 0;
}
题解：
终于会了反演= = 其实不就是一个公式吗= =....我们来膜拜PoPoQQQ （http://wenku.baidu.com/view/fbec9c63ba1aa8114431d9ac.html）
手推了半小时啊= =在推$f(i)$的时候竟然还不知道有约数和的公式这个玩意QAQ，感谢PoPoQQQ对我的教导，然后我根据这个玩意直接线性筛和mu一起预处理出来了...
首先我们变化一下公式，且先不考虑$<=a$这个条件，且下边均默认$n<=m$：
设
$$g(d) = \sum_{i=1}^{n} \sum_{j=1}^{m} [(i, j)=d]$$
$$f(d) = \sum_{i|d} i$$
那么本题可以变为枚举$gcd$，求
$$\sum_{i=1}^{n} g(i)f(i)$$
我们先来求$g(d)$（原来的博文是直接展开求的，现在我们用莫比乌斯反演来做...）
设
$$F(d)表示d|(i, j)的对数，而g(d)是表示(i, j)=d的对数$$
易得
$$ F(d)=
\begin{cases}
\lfloor \frac{n}{d} \rfloor \lfloor \frac{m}{d} \rfloor \\
\sum_{d|i} g(i) \\
\end{cases}
$$
根据莫比乌斯反演（两种形式，约数和 与 倍数和，这里是倍数和的形式，证明略），易得
$$ g(d) = \sum_{d|i} \mu (\frac{i}{d}) F(i) $$
$f(i)$的求法待会再说。那么现在我们展开之前得到的式子
$$\sum_{i=1}^{n} g(i)f(i) = \sum_{i=1}^{n} f(i) \sum_{i|T} \mu (\frac{T}{i}) F(T) = \sum_{i=1}^{n} f(i) \sum_{i|T} \mu (\frac{T}{i}) \lfloor \frac{n}{T} \rfloor \lfloor \frac{m}{T} \rfloor $$
继续换指标，我们来枚举$T$，最终可得
$$\sum_{T=1}^{n} \lfloor \frac{n}{T} \rfloor \lfloor \frac{m}{T} \rfloor \sum_{i|T} f(i) \mu (\frac{T}{i})$$
完美的公式= =
那么我们只需要求出$f(i)$然后按照倍数暴力更新即可，$O(nlnn)$（其实如果本题没有$<=a$这个条件，可以直接线性筛筛出来而不需要暴力枚举，详细请见我之前的博文http://www.cnblogs.com/iwtwiioi/p/4132095.html
那么考虑当有$<=a$这个条件了，发现其实就是$f(i)<=a$，那么发现其实只需要排序一下，然后用个bit动态维护前缀和即可，最终查询用分块。
复杂度$O(nlogn+qn^{0.5}log n)$
最后来说$f(i)$的线性求法...
由于
$$f(i) = \prod_{i} \sum_{j=0}^{a[i]} p_{i}^{j}，a[i]表示指数数目$$
我们考虑$f(kp_y)$，即线性筛里边的外层循环$k$和内层循环$p_y$
当$p_y | k$时
$$f(kp_y) = \prod_{i} \sum_{j=0}^{a[i]} p_{i}^{j}$$
我们在$p_y$那个指数提取一个出来，则那一部分的和变成$p_y \sum_{j=-1}^{a[y]-1} p_{y}^{j}$，将指数为$-1$的一起提出来，最终整理得
$$f(kp_y) = \prod_{i \neq y} \sum_{j=0}^{a[i]} p_{i}^{j} + p_y \prod_{i} \sum_{j=0}^{a[i]或当i=y时, a[y]-1} p_{i}^{j} = f(k去除所有p_y因子) + p_y f(k)$$
当$p_y \nmid k$时，好算多了...直接可推得
$$f(kp_y) = (1+p_y)f(k)$$
那么就大胆放到线性筛里去求吧！
