[bzoj4815][CQOI2017]小Q的表格

27 篇文章 0 订阅
24 篇文章 0 订阅

题目描述

小Q是个程序员。
作为一个年轻的程序员,小Q总是被老C欺负,老C经常把一些麻烦的任务交给小Q来处理。每当小Q不知道如何解决
时,就只好向你求助。为了完成任务,小Q需要列一个表格,表格有无穷多行,无穷多列,行和列都从1开始标号。
为了完成任务,表格里面每个格子都填了一个整数,为了方便描述,小Q把第a行第b列的整数记为f(a,b),为了完成
任务,这个表格要满足一些条件:(1)对任意的正整数a,b,都要满足f(a,b)=f(b,a);(2)对任意的正整数a,b,都要
满足b×f(a,a+b)=(a+b)*f(a,b)。为了完成任务,一开始表格里面的数很有规律,第a行第b列的数恰好等于a*b,
显然一开始是满足上述两个条件的。为了完成任务,小Q需要不断的修改表格里面的数,每当修改了一个格子的数
之后,为了让表格继续满足上述两个条件,小Q还需要把这次修改能够波及到的全部格子里都改为恰当的数。由于
某种神奇的力量驱使,已经确保了每一轮修改之后所有格子里的数仍然都是整数。为了完成任务,小Q还需要随时
获取前k行前k列这个有限区域内所有数的和是多少,答案可能比较大,只需要算出答案mod1,000,000,007之后的结
果。

性质

1、f(a,b)=f(b,a)
2、f(a,b)*(b-a)=f(a,b-a)*b
即f(a,b)/b=f(a,b-a)/(b-a)
考虑影响情况。
(a,b)<=>(b,a)
(a,b)<=>(a,b-a)
这是个辗转相除。
因此我们容易发现gcd相同的值会互相影响。
注意到对于(a,x),(a,b)和(y,b)互相的关系(gcd(a,x)=gcd(a,b)=gcd(y,b))。
f(a,x)/x=f(a,b)/b
f(a,b)/a=f(y,b)/y
容易发现对于gcd(x,y)=gcd(a,b)一定有f(x,y)/xy=f(a,b)/ab。
于是只需要保存a[i]=f(i,i)。
那么f(id,jd)=a[d]*ij其中(i,j)=1。
因此答案很容易计算,修改也很容易修改。
nd=1a[d]n/di=1n/dj=1[(i,j)=1]ij
后面部分可以设为g[n/d],考虑如何预处理它。

莫比乌斯意义

ni=1nj=1[(i,j)=1]ij
ni=1nj=1ijd|i,d|jμ(d)
nd=1μ(d)d2S(n/d)2
其中 S(n)=ni=1i
这个玩意可以n log n预处理,但还是不够快啊?
考虑作差,g[n]-g[n-1]是啥呢?
我们找到所有n/d!=(n-1)/d的d,那么d一定整除n。
d|nμ(d)d2[S(n/d)2S(n/d1)2]
d|nμ(d)d2[S(n/d)+S(n/d1)][S(n/d)S(n/d1)]
d|nμ(d)d2[2S(n/d)n/d]n/d
d|nμ(d)d2(n/d)3
n2d|nμ(d)n/d
n2ϕ(n)
这个为啥变成phi了?你考虑容斥意义理解即可。
于是线性筛phi,然后前缀和。

欧拉意义

也可以直接欧拉意义下理解。
考虑i>j,(i,j)=1。
我们知道当i>1时,所有和i互质的数两两配对,每两个的和都是i。
所以贡献和为 iiϕ(i)/2
注意到在j的时候也要算i,但是我们只考虑i>j,所以贡献乘2,就是
iiϕ(i)
i=1也恰好满足。

求答案

现在处理了g,可以每次分块计算。
我们需要支持修改和区间求和。
注意到修改m次,而求和有m sqrt n 次。
因此我们用分块维护,修改O(sqrt n),而询问O(1)。
即可平衡。

#include<cstdio>
#include<algorithm>
#define fo(i,a,b) for(i=a;i<=b;i++)
using namespace std;
typedef long long ll;
const int maxn=4000000+10,B=2000,mo=1000000007;
int sum[B+10],num[maxn],a[maxn],belong[maxn];
int pri[maxn],phi[maxn],g[maxn];
bool bz[maxn];
int i,j,k,d,y,l,t,n,m,p,tot,top,ans;
ll x;
int gcd(int a,int b){
    return b?gcd(b,a%b):a;
}
int getsum(int x){
    if (!x) return 0;
    int i,t=0,y=belong[x];
    t=sum[y-1];
    (t+=num[x])%=mo;
    return t;
}
int qsm(int x,int y){
    if (!y) return 1;
    int t=qsm(x,y/2);
    t=(ll)t*t%mo;
    if (y%2) t=(ll)t*x%mo;
    return t;
}
int main(){
    scanf("%d%d",&m,&n);
    fo(i,1,n){
        belong[i]=(i-1)/B+1;
        a[i]=(ll)i*i%mo;
        if (belong[i]==belong[i-1]){
            num[i]=(num[i-1]+a[i])%mo;
            (sum[belong[i]]+=a[i])%=mo;
        }
        else{
            num[i]=a[i];
            sum[belong[i]]=(sum[belong[i]-1]+a[i])%mo;
        }
    }
    phi[1]=1;
    fo(i,2,n){
        if (!bz[i]){
            pri[++top]=i;
            phi[i]=i-1;
        }
        fo(j,1,top){
            if ((ll)i*pri[j]>n) break;
            bz[i*pri[j]]=1;
            if (i%pri[j]==0){
                phi[i*pri[j]]=phi[i]*pri[j];
                break;
            }
            phi[i*pri[j]]=phi[i]*(pri[j]-1);
        }
    }
    fo(i,1,n) g[i]=(g[i-1]+(ll)i*i%mo*phi[i]%mo)%mo;
    while (m--){
        scanf("%d%d%lld%d",&j,&k,&x,&p);
        t=x%mo;
        d=gcd(j,k);
        t=(ll)t*qsm(j/d,mo-2)%mo*qsm(k/d,mo-2)%mo;
        a[d]=t;
        y=belong[d];
        l=0;
        fo(i,(y-1)*B+1,min(y*B,n)){
            num[i]=(l+a[i])%mo;
            l=num[i];
        }
        (l+=sum[y-1])%=mo;
        l=(l-sum[y])%mo;
        fo(i,y,belong[n]) (sum[i]+=l)%=mo;
        ans=0;
        i=1;
        while (i<=p){
            j=p/(p/i);
            t=(getsum(j)-getsum(i-1))%mo;
            (ans+=(ll)t*g[p/i]%mo)%=mo;
            i=j+1;
        }
        (ans+=mo)%=mo;
        printf("%d\n",ans);
    }
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值