http://www.elijahqi.win/archives/3477
求∑∑((n mod i)*(m mod j))其中1<=i<=n,1<=j<=m,i≠j。
这份代码调试了很久..菜
把式子随便推导一下展开即可 注意取模 6的逆元直接exgcd算完
公式:
∑ni=1∑mj=1(nmodi)×(mmodj)−∑min(n,m)i=1(nmodi)×(mmodj) ∑ i = 1 n ∑ j = 1 m ( n mod i ) × ( m mod j ) − ∑ i = 1 m i n ( n , m ) ( n mod i ) × ( m mod j ) ∑ni=1∑mj=1(n−i×⌊ni⌋)×(m−j×⌊mj⌋)−∑min(n,m)i=1(n−i×⌊ni⌋)×(m−j×⌊mj⌋) ∑ i = 1 n ∑ j = 1 m ( n − i × ⌊ n i ⌋ ) × ( m − j × ⌊ m j ⌋ ) − ∑ i = 1 m i n ( n , m ) ( n − i × ⌊ n i ⌋ ) × ( m − j × ⌊ m j ⌋ ) ∑ni=1(n−i×⌊ni⌋)∑mj=1(m−j×⌊mj⌋)−∑min(n,m)i=1(n−i×⌊ni⌋)×(m−j×⌊mj⌋)) ∑ i = 1 n ( n − i × ⌊ n i ⌋ ) ∑ j = 1 m ( m − j × ⌊ m j ⌋ ) − ∑ i = 1 m i n ( n , m ) ( n − i × ⌊ n i ⌋ ) × ( m − j × ⌊ m j ⌋ ) )
#include<cstdio>
#include<algorithm>
#define ll long long
using namespace std;
const int mod=19940417;
const int inv=3323403;
inline int inc(int x,int v){return x+v>=mod?x+v-mod:x+v;}
inline int dec(int x,int v){return x-v<0?x-v+mod:x-v;}
inline int gao(int n,int m){
int tmp=0,last;
for (int i=1;i<=n;i=last+1){
last=m/(m/i);if (last>n) last=n;
tmp=inc(tmp,(ll)(m/i)*((ll)(last-i+1)*(last+i)/2%mod)%mod);
}return tmp;
}
inline int calc(int x){return (ll)x*(x+1)%mod*(2*x+1)%mod*inv%mod;}
int n,m,last;
int main(){
freopen("bzoj2956.in","r",stdin);
scanf("%d%d",&n,&m);int ans=0,ans1=0,ans2=0;
ans1=inc(ans1,(ll)n*n%mod);ans2=inc(ans2,(ll)m*m%mod);
ans1=dec(ans1,gao(n,n));ans2=dec(ans2,gao(m,m));
if (n>m) swap(n,m);ans=inc(ans,(ll)n*m%mod*n%mod);
for (int i=1;i<=n;i=last+1){
last=min(n/(n/i),m/(m/i));
ans=inc(ans,(ll)(m/i)*(n/i)%mod*dec(calc(last),calc(i-1))%mod);
}
// printf("%d\n",m*gao(n,n));
// printf("%d\n",n*gao(n,m));
ans=dec(ans,(ll)m*gao(n,n)%mod);ans=dec(ans,(ll)n*gao(n,m)%mod);
ans1=dec((ll)ans1*ans2%mod,ans);printf("%d\n",ans1);
return 0;
}