【数论分块】[BZOJ2956、LuoguP2260] 模积和

十年OI一场空，忘记取模见祖宗

题目：

求$$\sum_{i=1}^{n}\sum_{j=1}^{m} (n \bmod i)(m \bmod i)$$

(其中i,j不相等)

暴力拆式子：

　$$ANS=\sum_{i=1}^{n}\sum_{j=1}^{m} (n- \left \lfloor \frac{n}{i} \right \rfloor*i)(m- \left \lfloor \frac{m}{i} \right \rfloor*i)-\sum_{i=1}^{min(n,m)} (n- \left \lfloor \frac{n}{i} \right \rfloor *i)(m- \left \lfloor \frac{m}{i} \right \rfloor *i)$$

令$f(n)=\sum_{i=1}^{n} (n- \left \lfloor \frac{n}{i} \right \rfloor *i)$

令$g(n)=\sum_{i=1}^{n}(n- \left \lfloor \frac{n}{i} \right \rfloor *i)(m- \left \lfloor \frac{m}{i} \right \rfloor *i)$

不妨设n<=m

则有

$$ANS=f(n)*f(m)-g(n)$$

其中$$g(n)=\sum_{i=1}^{n} n*m-n*\sum_{i=1}^{n} \left \lfloor \frac{m}{i} \right \rfloor *i-m*\sum_{i=1}^{n} \left \lfloor \frac{n}{i} \right \rfloor *i+\sum_{i=1}^{n} \left \lfloor \frac{n}{i} \right \rfloor* \left \lfloor \frac{m}{i} \right \rfloor *i^2$$

且易有$$\sum_{i=1}^{n} i^2=\frac{n*(n+1)*(2*n+1)}{6}$$

预处理6在模19940417意义下的逆元（我用了exgcd）

然后用数论分块把上面一堆东西算一下即可

 #include<bits/stdc++.h>

 #define int long long

 #define writeln(x)  write(x),puts("")

 #define writep(x)   write(x),putchar(' ')

 using namespace std;

 inline int read(){

     int ans=,f=;char chr=getchar();

     while(!isdigit(chr)){if(chr=='-') f=-;chr=getchar();}

     while(isdigit(chr)){ans=(ans<<)+(ans<<)+chr-;chr=getchar();}

     return ans*f;

 }void write(int x){

     if(x<) putchar('-'),x=-x;

     if(x>) write(x/);

     putchar(x%+'');

 }const int mod = ;

 int n,m,k;

 inline void Add(int &x,int y){x+=y;x%=mod;}

 void exgcd(int a,int b,int &x,int &y){

     if(b==)return x=,y=,void();

     exgcd(b,a%b,x,y);

     int t=x;x=y,y=t-a/b*y;

 }int inv(int x){

     int xx,y;

     exgcd(,mod,xx,y);

     xx=(xx%mod+mod)%mod;

     return xx;

 }const int inv6=inv();

 int sum(int x){return (x)*(x+)%mod*(*x%mod+)%mod*inv6%mod;}

 int query1(int l,int r){return ((sum(r)-sum(l-))%mod+mod)%mod;}

 int query2(int l,int r){int ans=(r-l+)*(l+r)/;return ans%mod;}

 int calc1(int n){

     int ans=;

     for(int i=,j,t;i<=n;i=j+){

         j=n/(n/i);

         t=n/i*(i+j)*(j-i+)/;

         t%=mod;

         Add(ans,t);

     }ans=n*n%mod-ans;

     ans=(ans%mod+mod)%mod;

     return ans;

 }int calc2(int k){

     int ans=;

     for(int i=,j,t;i<=n;i=j+){

         j=min(n/(n/i),m/(m/i));

         int s1=n*(m/i)%mod*query2(i,j)%mod;

         int s2=m*(n/i)%mod*query2(i,j)%mod;

         int s3=(n/i)*(m/i)%mod*query1(i,j)%mod;

         Add(s1,s2);

         Add(ans,s1);

         ans-=s3;

         ans=((ans)%mod+mod)%mod;

     }return ans;

 }

 signed main(){

     n=read(),m=read();

     if(n>m)swap(n,m);

     int ans=calc1(n)*calc1(m)%mod;

     ans-=n*m%mod*n%mod;

     ans=(ans%mod+mod)%mod;

     ans+=calc2(n);

     ans=(ans%mod+mod)%mod;

     cout<<ans<<endl;

     return ;

 }

巴特西

【数论分块】[BZOJ2956、LuoguP2260] 模积和

最新文章

热门文章