[HDU4947]GCD Array(莫比乌斯反演)
题面
http://acm.hdu.edu.cn/showproblem.php?pid=4947
前置知识
- 线性筛积性函数：https://www.cnblogs.com/zhoushuyu/p/8275530.html
- 莫比乌斯反演：《具体数学：计算机科学基础》4.9节
- 莫比乌斯反演的基础应用： https://www.cnblogs.com/xh092113/p/12268356.html
- 数论分块： https://www.cnblogs.com/xh092113/p/12269753.html
题解
我们可以维护数组f，使得f时刻等于\(a{\times}{\mu}\)。
考虑每次操作1 n d v对于数组a的改变量\({\Delta}a(x)\)，发现
\[{\Delta}a(x)=v*[(x,n)=d]
\]
\[=v[d|x][d|n][({\frac{x}{d}},{\frac{n}{d}})=1]
\]
如果\(d{\nmid}n\)，那么这个操作可以直接忽略。所以不妨设d|n，
\[{\Delta}a(x)=v[d|x][({\frac{x}{d}},{\frac{n}{d}})=1]
\]
利用消gcd的技巧，反演得
\[{\Delta}a(x)=v[d|x]{\sum\limits_{p|{\frac{x}{d}},p|{\frac{n}{d}}}}{\mu}(p)
\]
\[={\sum\limits_{p|{\frac{n}{d}}}}v{\mu}(p)[d|x][p|{\frac{x}{d}}]
\]
\[={\sum\limits_{p|{\frac{n}{d}}}}v{\mu}(p)[pd|x]
\]
此时的\({\Delta}f\)也应该等于\({\Delta}a{\times}{\mu}\)。
\[{\Delta}f(x)={\sum\limits_{y|x}}{\Delta}a({\frac{x}{y}}){\mu(y)}
\]
\[={\sum\limits_{y|x}}{\mu}(y){\sum\limits_{p|{\frac{n}{d}}}}v{\mu}(p)[pd|{\frac{x}{y}}]
\]
\[={\sum\limits_{p|{\frac{n}{d}}}}v{\mu(p)}{\sum\limits_{y|x}}{\mu(y)}[pd|{\frac{x}{y}}]
\]
\[={\sum\limits_{p|{\frac{n}{d}}}} v{\mu(p)} [pd|x] {\sum\limits_{y|{\frac{x}{pd}}}} {\mu}(y)
\]
\[={\sum\limits_{p|{\frac{n}{d}}}} v{\mu(p)} [pd|x] [{\frac{x}{pd}}=1]
\]
\[={\sum\limits_{p|{\frac{n}{d}}}}v{\mu(p)}[x=pd]
\]
由此可见，对于1 n d v的修改操作，只需要遍历\({\frac{n}{d}}\)的所有因子p，然后将\(f(pd)+=v{\mu(p)}\)即可。
对于查询操作，同样将结果用f表示：
\[{\sum\limits_{i=1}^{x}}a(i)={\sum\limits_{i=1}^{x}}{\sum\limits_{p|i}}f(p)
\]
\[={\sum\limits_{p}}{\sum\limits_{i=1}^{x}}[p|i]f(p)
\]
\[={\sum\limits_{p}}{\lfloor}{\frac{x}{p}}{\rfloor}f(p)
\]
将f的前缀和使用树状数组维护，数论分块即可优化到每次询问\(O({\sqrt{x}}log x)\)。
但是，注意到这里的区间查询是O(logx)而非普通数论分块的O(1)，因此此处可以进一步优化：
设置阈值T，当查询的x小于T时使用暴力求和，否则使用数论分块。
则暴力求和的时间复杂度为O(T)，数论分块的时间复杂度为\(O({\frac{x}{T}}log x)\)。不难发现置\(T={\sqrt{xlog x}}\)时总时间复杂度最低。至此，每次询问的时间复杂度降为\(O(\sqrt{x log x})\)。
分析总时间，各预处理均在nlogn级别及以下，（n为数组长度），每次修改是\(O({\sqrt{n}} log n)\)，查询是\(O(\sqrt{n log n})\)。故总时间复杂度应为\(O(CQ{\sqrt{n}} log n)\)，（C为数据组数)，且此处的根号来源于因子个数，其常数极小。
代码
#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define rg register
#define N 50000
#define W 200000
inline ll read(){
ll s = 0,ww = 1;
char ch = getchar();
while(ch < '0' || ch > '9'){if(ch == '-')ww = -1;ch = getchar();}
while('0' <= ch && ch <= '9'){s = 10 * s + ch - '0';ch = getchar();}
return s * ww;
}
inline void write(ll x){
if(x < 0)putchar('-'),x = -x;
if(x > 9)write(x / 10);
putchar('0' + x % 10);
}
ll n,Q;
ll a[N+5],f[N+5],lg[W+5];
vector<ll>fac[W+5];
inline void prepro(){ //预处理因子和log
for(rg ll i = 1;i <= W;i++)
for(rg ll j = i;j <= W;j += i)fac[j].push_back(i);
lg[1] = 0;
for(rg ll i = 2;i <= W;i++)lg[i] = lg[i>>1] + 1;
}
ll pn;
bool isp[W+5];
ll pri[40000+5],mu[W+5];
inline void Eular(){ //线性筛
mu[1] = 1;
for(rg ll i = 2;i <= W;i++)isp[i] = 1;
for(rg ll i = 2;i <= W;i++){
if(isp[i])pri[++pn] = i,mu[i] = -1;
for(rg ll j = 1;i * pri[j] <= W;j++){
isp[i*pri[j]] = 0;
if(i % pri[j])mu[i*pri[j]] = -mu[i];
else{
mu[i*pri[j]] = 0;
break;
}
}
}
}
inline ll lowbit(ll x){
return x & -x;
}
struct BIT{
ll b[N+5];
inline void reset(){
memset(b,0,sizeof(b));
}
inline void ud(ll x,ll dx){
for(rg ll i = x;i <= n;i += lowbit(i))b[i] += dx;
}
inline ll query(ll x){
ll rt = 0;
for(rg ll i = x;i;i -= lowbit(i))rt += b[i];
return rt;
}
inline ll sum(ll l,ll r){
return query(r) - query(l - 1);
}
}B;
int main(){
prepro();
Eular();
n = read(),Q = read();
ll t = 0;
while(n || Q){
printf("Case #%lld:\n",++t);
B.reset();
memset(a,0,sizeof(a));
memset(f,0,sizeof(f));
while(Q--){
ll opt = read();
if(opt == 1){
ll _n = read(),d = read(),v = read();
if(_n % d)continue;
for(rg ll i = 0;i < fac[_n/d].size();i++){
ll p = fac[_n/d][i];
if(p * d > n)break;
f[p*d] += mu[p] * v;
B.ud(p * d,mu[p] * v);
}
}
else{
ll x = read();
ll T = (ll)round(sqrt(x*lg[x]));
if(T > x)T = x;
ll ans = 0;
for(rg ll i = 1;i <= T;i++)ans += f[i] * (x / i); //sqrt(xlogx)以下暴力求和
ll L,R = T;
while(R != x){ //sqrt(xlogx)以上数论分块
L = R + 1;
R = x / (x / L);
ans += B.sum(L,R) * (x / L);
}
write(ans),putchar('\n');
}
}
n = read(),Q = read();
}
return 0;
}
