【jzoj3117】【WinterCamp 2013】【模积和】【分块】【扩展欧几里得】

题目大意

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

解题思路

观察可知ans=∑∑((n mod i)(m mod j))1<=i<=n,1<=j<=m,i≠j。

=>ans=∑(n mod i)∑(m mod j)-∑((n mod i)(m mod i))

f(n)=∑(n mod i)

g(n,m)=∑(n mod i)(m mod i)

ans=f(n)f(m)-g(n,m)

f(n)=∑(n−⌊n/i⌋i)

对于⌊n/i⌋我们可以使用分块，对于[i,j]满足n/i==n/j可以一起处理，对于n可以用乘法求，对于*i可以用等差数列求和。

g(n,m)=∑(nm+⌊n/i⌋⌊m/i⌋)i2−⌊n/i⌋mi−⌊m/i⌋ni)

对于第一项可以用乘法解决，对于第二项可以分块后用完全平方和求和公式求，第三四项用等差数列求和。

∑(1^2+2^2+……+n^2)==n(n+1)(2n+1)/6详细证明不再给出

观察可知模数不是质数，使用扩展欧几里得求出逆元。

aa’==1(%mo)

=>aa’+moy==1

由于a，mo互质，a’可以用扩展欧几里得求出。

code

#include<cmath>#include<cstdio>#include<cstring>#include<algorithm>#define LF double#define LL long long#define min(a,b) ((a<b)?a:b)#define max(a,b) ((a>b)?a:b)#define fo(i,j,k) for(int i=j;i<=k;i++)#define fd(i,j,k) for(int i=j;i>=k;i--)using namespace std;int const mo=19940417;int n,m,ni;int calc(int n){    int ans=1ll*n*n%mo;    for(int i=1,j=1;i<=n;i=j+1){        j=n/(n/i);        ans=(ans-1ll*(i+j)*(j-i+1)/2%mo*(n/i)%mo)%mo;    }    return (ans+mo)%mo;}int calc3(int n){    return 1ll*n*(n+1)%mo*(2*n+1)%mo*ni%mo;}int calc2(int n,int m){    int mi=min(n,m),ans=1ll*n*m%mo*mi%mo;    for(int i=1,j=1;i<=mi;i=j+1){        j=min(n/(n/i),m/(m/i));        ans=(ans+1ll*(calc3(j)-calc3(i-1))%mo*(n/i)%mo*(m/i)%mo                      -1ll*(i+j)*(j-i+1)/2%mo*(n/i)%mo*m%mo                      -1ll*(i+j)*(j-i+1)/2%mo*(m/i)%mo*n%mo)%mo;    }    return (ans+mo)%mo;}void exgcd(int a,int &x,int b,int &y){    if(!b){x=1;y=0;return;}    exgcd(b,y,a%b,x);    y-=a/b*x;}int main(){    freopen("d.in","r",stdin);    freopen("d.out","w",stdout);    scanf("%d%d",&n,&m);    int y;exgcd(6,ni,mo,y);    printf("%d",(1ll*calc(n)*calc(m)%mo-calc2(n,m)+mo)%mo);    return 0;}