zoukankan      html  css  js  c++  java
  • 【BZOJ3529】[SDOI2014] 数表(莫比乌斯反演)

    点此看题面

    大致题意: 规定一个(n*m)数表中每个数为(sum_{d|i,d|j}d),求数表中不大于(a)的数之和。

    不考虑限制

    我们先不考虑限制,来推一波式子。

    首先,易知数表中第(i)行第(j)列的数应该是(sigma(gcd(i,j)))

    则和就为:

    [sum_{d=1}^{min(n,m)}sigma(d)sum_{i=1}^{lfloorfrac nd floor}sum_{j=1}^{lfloorfrac md floor}[gcd(i,j)=1] ]

    ([gcd(i,j)=1])可以化成(sum_{p|gcd(i,j)}mu(p)),若枚举(p),就得到:

    [sum_{d=1}^{min(n,m)}sigma(d)sum_{p=1}^{lfloorfrac{min(n,m)}d floor}mu(p)lfloorfrac n{dp} floorlfloorfrac m{dp} floor ]

    (g=dp),调整枚举顺序得到:

    [sum_{g=1}^{min(n,m)}lfloorfrac n{dp} floorlfloorfrac m{dp} floorsum_{d|g}sigma(d)mu(frac gd) ]

    离线处理限制

    考虑上面的式子只有当(sigma(d)le a)时才会被计算答案。

    则我们考虑设(T(g)=sum_{d|g}sigma(d)mu(frac gd)),一开始全为(0)

    然后我们按照(a)从小到大枚举询问,每次将(sigma(d)le a)(d)(10^5)范围内的倍数所对应的(T(g))全都加上(sigma(d)mu(frac gd))

    但注意到询问时使用除法分块需要求一段区间的(T)值和,则我们用树状数组维护就可以了。

    关于取模的细节

    注意,这里的取模是向(2^{31})取模,则我们可以考虑先开(unsigned int)计算答案,最后再将其向(2^{31})取模

    代码

    #include<bits/stdc++.h>
    #define Tp template<typename Ty>
    #define Ts template<typename Ty,typename... Ar>
    #define Reg register
    #define RI Reg int
    #define Con const
    #define CI Con int&
    #define I inline
    #define W while
    #define N 100000
    #define Q 20000
    #define MxS 500000
    #define UI unsigned
    #define RU Reg unsigned
    #define CU Con unsigned&
    #define LL long long
    #define Gmax(x,y) (x<(y)&&(x=(y)))
    #define min(x,y) ((x)<(y)?(x):(y))
    #define max(x,y) ((x)>(y)?(x):(y))
    #define pb push_back
    #define IT vector<int>::iterator
    using namespace std;
    int Qt,Qans[Q+5];vector<int> s[MxS+5];
    struct Query//询问
    {
    	int x,y,v,pos;I Query(CI a=0,CI b=0,CI z=0,CI p=0):x(a),y(b),v(z),pos(p){}
    	I bool operator < (Con Query& o) Con {return v<o.v;}
    }q[Q+5];
    class FastIO
    {
    	private:
    		#define FS 100000
    		#define tc() (A==B&&(B=(A=FI)+fread(FI,1,FS,stdin),A==B)?EOF:*A++)
    		#define pc(c) (C==E&&(clear(),0),*C++=c)
    		#define tn (x<<3)+(x<<1)
    		#define D isdigit(c=tc())
    		int T;char c,*A,*B,*C,*E,FI[FS],FO[FS],S[FS];
    	public:
    		I FastIO() {A=B=FI,C=FO,E=FO+FS;}
    		Tp I void read(Ty& x) {x=0;W(!D);W(x=tn+(c&15),D);}
    		Ts I void read(Ty& x,Ar&... y) {read(x),read(y...);}
    		Tp I void write(Ty x) {W(S[++T]=x%10+48,x/=10);W(T) pc(S[T--]);}
    		Tp I void writeln(Con Ty& x) {write(x),pc('
    ');}
    		I void clear() {fwrite(FO,1,C-FO,stdout),C=FO;}
    }F;
    class LinearSiever//线性筛
    {
    	private:
    		int Pt,P[N+5],Mn[N+5];
    		I LL Qpow(LL x,LL y) {LL t=1;W(y) y&1&&(t*=x),x*=x,y>>=1;return t;}//快速幂
    	public:
    		int MxSigma,sigma[N+5],mu[N+5];
    		I void Sieve(CI S)
    		{
    			RI i,j,x,t;for(mu[1]=1,i=2;i<=S;++i)//筛mu,筛最小质因数用于求sigma
    			{
    				!Mn[i]&&(mu[P[++Pt]=i]=-1,Mn[i]=i);
    				for(j=1;j<=Pt&&1LL*i*P[j]<=S;++j)
    					if(Mn[i*P[j]]=P[j],i%P[j]) mu[i*P[j]]=-mu[i];else break;
    			}
    			for(sigma[1]=1,i=2;i<=S;++i)//求sigma
    			{
    				x=i,t=0;W(!(x%Mn[i])) x/=Mn[i],++t;
    				sigma[i]=sigma[x]*((Qpow(Mn[i],t+1)-1)/(Mn[i]-1)),Gmax(MxSigma,sigma[i]);
    			}
    		}
    }L;
    class TreeArray//树状数组
    {
    	private:
    		#define lowbit(x) (x&-x)
    		UI v[MxS+5];
    		I UI QS(RI x) {RU t=0;W(x) t+=v[x],x-=lowbit(x);return t;}//询问前缀
    	public:
    		I void Add(RI x,CI y) {W(x<=L.MxSigma) v[x]+=y,x+=lowbit(x);}//单点修改
    		I UI Qry(CI l,CI r) {return QS(r)-QS(l-1);}//区间查询
    }T;
    I void Upt(CI x,CI v) {for(RI i=1;1LL*x*i<=N;++i) T.Add(x*i,L.sigma[x]*L.mu[i]);}//更新一个数倍数的值
    int main()
    {
    	RI i,p=1,t,x,y,v,l,r;UI ans;for(L.Sieve(N),i=1;i<=N;++i) s[L.sigma[i]].pb(i);//用桶对sigma值进行排序
    	for(F.read(Qt),i=1;i<=Qt;++i) F.read(x,y,v),q[i]=Query(min(x,y),max(x,y),v,i);//读入询问
    	for(sort(q+1,q+Qt+1),i=1;i<=Qt;++i)//对询问按a从小到大排序
    	{
    		W(p<=q[i].v) {for(IT it=s[p].begin();it!=s[p].end();++it) Upt(*it,p);++p;}//更新sigma(d)≤a的d的倍数的T值
    		for(ans=0,t=min(q[i].x,q[i].y),l=1;l<=t;l=r+1)//除法分块
    			r=min(q[i].x/(q[i].x/l),q[i].y/(q[i].y/l)),ans+=T.Qry(l,r)*(q[i].x/l)*(q[i].y/l);
    		Qans[q[i].pos]=ans%(1LL<<31);//存储答案并取模
    	}
    	for(i=1;i<=Qt;++i) F.writeln(Qans[i]);return F.clear(),0;//输出答案
    }
    
  • 相关阅读:
    我用自己做的图书比价搜索买了一本书
    2.17
    最近的工作
    FireBug的Bug
    2.18
    tecent面试题解答
    .net杂记
    python的round测试
    最近在网上买书的体会
    关于迅雷评论的一个改造html css
  • 原文地址:https://www.cnblogs.com/chenxiaoran666/p/BZOJ3529.html
Copyright © 2011-2022 走看看