題面
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;
}