Bzoj3853 GCD Array

Time Limit: 6 Sec  Memory Limit: 64 MB
Submit: 328  Solved: 125

Description

Teacher Mai finds that many problems about arithmetic function can be reduced to the following problem:

Maintain an array a with index from 1 to l. There are two kinds of operations:

  1. Add v to ax for every x that gcd(x,n)=d.

  2. Query Sigma(Xi) (i<=1<=X)

Input

There are multiple test cases, terminated by a line "0 0".

For each test case, the first line contains two integers l,Q(1<=l,Q<=5*10^4), indicating the length of the array and the number of the operations.

In following Q lines, each line indicates an operation, and the format is "1 n d v" or "2 x" (1<=n,d,v<=2*10^5,1<=x<=l).

Output

For each case, output "Case #k:" first, where k is the case number counting from 1.

Then output the answer to each query.

 

Sample Input

6 4
1 4 1 2
2 5
1 3 3 3
2 3
0 0

Sample Output

Case #1:
6
7

HINT

Source

By 镇海中学

数学问题 莫比乌斯反演

我们注意到每次修改的位置是不连续的，很难用数据结构统一维护。但是这个数据范围看上去必须用数据结构维护，那么我们就得想办法让每次修改/查询的部分连续起来。

方法一：分块。

　　将修改的数量分块，如果$ gcd(x,n)=d $的个数大于 $ sqrt(i) $个，暴力修改，如果小于$sqrt(i)$，分块修改，有一部分带修改点构成了有规律的区间。

　　常数优秀的话应该可以过（口胡，没尝试过）

方法二：莫比乌斯反演

　　我们把题目中各个关键条件都写在草稿纸上，这样有很大几率会看到类似$ v [gcd(x,n)=d]$ 的式子

　　看到这个形式，我的内心充满波动，甚至想要反演！

　　试着推一下式子：
$$ v[gcd(x,n)=d] $$
$$ =v[gcd(x/d,n/d)=1] $$
$$ =v * sum_{k|frac{x}{d} , k|frac{n}{d}} mu(k)$$
把
    $[k|frac{x}{d}] mu(k)$
变形成
    $[kd|x] mu(k)$
那么我们就可以O($sqrt n$)地枚举 $frac{n}{d}$的因数k，再单点修改每一个 $kd|x$ 的位置。

(x是不确定的位置，修改了每一个kd，就等于修改了所有受影响的x)
嗯？我们好像并没有解决单点修改的问题？
但我们可以区间查询了！
设原来要维护的数组是a[]，现在维护的数组是c[]，
那么有
$$a[i]=sum_{j|i}c[j]$$
于是
$$ans=sum_{i=i}^{x} c[i] =sum_{i=i}^{x} sum_{j|i}c[j] $$
$$ans=sum_{j=1}^{x}c[j]* lfloor frac{x}{j} floor$$
可以应用分块加速的套路。
用数据结构维护c[j]的前缀和即可。

一个脑洞：

　　如果每次给定k和b，修改满足"kx+b<=n"的"kx+b"位置，要怎么玩？(k,b,x都是整数）

　　“有问必答不知道”

另一个脑洞：

　　如果每次修改$ lcm(x,n)=d $ 的位置，要怎么玩？(x,n,d都是整数）

　　我们试着推一下：

　　$$ v[lcm(x,n)=d] $$
　　$$ v[gcd(x,n)=frac{x*n}{d} ]$$
　　$$ v[gcd(d/n,d/x)=1]$$
　　$$ =v * sum_{k|frac{d}{x} , k|frac{d}{n}} mu(k)$$

　　我们发现这次x到了分母的位置，且d每次修改的时候都会变，维护起来就很麻烦……
　　嗯……我编不下去了！（逃）

原题代码：

#include<iostream>
#include<algorithm>
#include<cstring>
#include<cstdio>
#include<cmath>
#include<vector>
#define LL long long
using namespace std;
const int mxn=200010;
int read(){
    int x=0,f=1;char ch=getchar();
    while(ch<'0' || ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0' && ch<='9'){x=x*10+ch-'0';ch=getchar();}
    return x*f;
}
int n,Q;
LL t[mxn];
void add(int x,int v){
    while(x<=n){t[x]+=v;x+=x&-x;}return;
}
LL smm(int x){
    LL res=0;
    while(x){res+=t[x];x-=x&-x;}
    return res;
}
//
int pri[mxn],mu[mxn],cnt=0;
bool vis[mxn];
vector<int>ve[mxn];
void init(){
    mu[1]=1;
    for(int i=2;i<mxn;i++){
        if(!vis[i]){
            pri[++cnt]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=cnt && pri[j]*i<mxn;j++){
            vis[pri[j]*i]=1;
            if(i%pri[j]==0){mu[i*pri[j]]=0;break;}
            mu[pri[j]*i]=-mu[i];
        }
    }
    for(int i=1;i<mxn;i++)
        for(int j=i;j<mxn;j+=i)
            ve[j].push_back(i);
    return;
}
LL calc(int x){
    LL res=0;
    for(int i=1,pos;i<=x;i=pos+1){
        int y=x/i; pos=x/y;
        res=res+(smm(pos)-smm(i-1))*(LL)y;
    }
    return res;
}
int main(){
    int i,j;
    init();
    int cas=0;
    while(scanf("%d%d",&n,&Q)!=EOF){
        if(!n && !Q)break;
        printf("Case #%d:
",++cas);
        memset(t,0,sizeof t);
        int op,x,v,d;
        while(Q--){
            op=read();
            if(op==1){
                x=read();d=read();v=read();
                if(x%d==0){
                    int st=x/d;
                    for(j=0;j<ve[st].size();j++){
                        int dx=ve[st][j];
                        add(dx*d,v*mu[dx]);
                    }
                }
            }
            else{
                x=read();
                LL ans=calc(x);
                printf("%lld
",ans);
            }
        }
    }
    return 0;
}