BZOJ 2956: 模积和

题目

2956: 模积和

Time Limit: 10 Sec  Memory Limit: 128 MB
Submit: 554  Solved: 257
[Submit][Status]

Description

　求∑∑((n mod i)*(m mod j))其中1<=i<=n,1<=j<=m,i≠j。

　　

Input

第一行两个数n,m。

Output

　　一个整数表示答案mod 19940417的值

Sample Input

3 4

Sample Output

1

样例说明
　　答案为(3 mod 1)*(4 mod 2)+(3 mod 1) * (4 mod 3)+(3 mod 1) * (4 mod 4) + (3 mod 2) * (4 mod 1) + (3 mod 2) * (4 mod 3) + (3 mod 2) * (4 mod 4) + (3 mod 3) * (4 mod 1) + (3 mod 3) * (4 mod 2) + (3 mod 3) * (4 mod 4) = 1

数据规模和约定
　　对于100%的数据n,m<=10^9。

题解

这道题目我调了好久，有两点错：一是乱用除法QAQ，要乘逆元才对！二是错误的估计了结果的大小，少用了几个%。虽然取模的数很小，但是两个乘在一起照样让你爆精度！每步都要取模！

这题的思路是分别计算∑∑((n mod i)*(m mod j))和∑∑((n mod i)*(m mod j))[i==j]的值，并且作差。

最后化简得式子是

∑∑((n mod i) * (m mod j)) 1<=i<=n, 1<=j<=m, i≠j=

∑(n mod i) * ∑(m mod i) - ∑((n mod i) * (m mod i))=

∑(n-[n/i]*i) * ∑(m-[m/i]*i) - ∑(nm-([n/i]+[m/i])i+[n/i][m/i]*i*i)【右式中[]表示下取整】

代码

/*Author:WNJXYK*/
#include<cstdio>
#include<iostream>
using namespace std;
const long long M=19940417;

long long n,m;


inline int remin(int a,int b){
    if (a<b) return a;
    return b;
}

inline long long sum(long long n){  
    return n*(n+1)%M*(2*n+1)%M*3323403%M;  
} 

inline long long getDist(int x,int y){
    return (x+y)%M*(y-x+1)%M*9970209%M;
} 

inline long long getOneAns(long long x){
    long long ans=((long long)x*(long long)x)%M;
    long long pos;
    for (long long i=1;i<=x;i=pos+1){
        pos=remin(x/(x/i),x);
        //cout<<"+"<<(long long)x*(long long)(pos-i+1)%M<<endl;
        //cout<<"-"<<((long long)((long long)(pos)*(long long)(pos+1)/(long long)2-(long long)(i-1)*(long long)(i)/(long long)2)*(long long)(x/i))%M<<endl;
        //cout<<"x"<<(long long)(pos)*(long long)(pos+1)/(long long)2-(long long)(i-1)*(long long)(i)/(long long)2<<endl;
        //ans=ans-((long long)getDist(i,pos)%M*(long long)(x/i))%M;
        ans-=(x/i)*getDist(i,pos)%M;
        while(ans<0) ans+=M;
    }
    //cout<<endl;
    while(ans<0) ans+=M;
    return ans;
} 

inline long long getTwoAns(long long x1,long long x2){
    long long ans=0;
    long long pos;
    long long x=remin(x1,x2);
    ans=x%M*x1%M*x2%M;
    for (long long i=1;i<=x;i=pos+1){
        pos=remin(x1/(x1/i),x2/(x2/i));
        //cout<<"POS"<<i<<"~"<<pos<<endl;
        //cout<<"+"<<(long long)(pos-i+1)*(long long)x1*(long long)x2<<endl;
        //cout<<"-"<<(long long)x1*(long long)(x2/i)*getDist(i,pos)<<endl;
        ///cout<<"-"<<(long long)x2*(long long)(x1/i)*getDist(i,pos)<<endl;
        //cout<<"+"<<((long long)((long long)x1/(long long)i)*(long long)((long long)x2/(long long)i))*getSDist(i,pos)<<endl;
        //cout<<"————————————"<<endl;
        ans=ans-(long long)x1%M*(long long)(x2/i)%M*getDist(i,pos)%M-(long long)x2%M*(long long)(x1/i)%M*getDist(i,pos)%M+((long long)((long long)x1/(long long)i)%M*(long long)((long long)x2/(long long)i))%M*(sum(pos)-sum(i-1))%M;
        while(ans<0) ans+=M;
    }
    while(ans<0) ans+=M;
    return ans;
}

int main(){
    scanf("%lld%lld",&n,&m);
    long long Ans=((long long)getOneAns(n)*(long long)getOneAns(m))%M;
    Ans-=getTwoAns(n,m);
    //cout<<getOneAns(n)<<endl;
    //cout<<getOneAns(m)<<endl;
    //cout<<getTwoAns(n,m)<<endl;
    while(Ans<0) Ans+=M;
    printf("%lld
",Ans%M);
    return 0;
}