P2260 [清华集训2012]模积和

由

(amod b=a-b(a/b))

所以

(sumlimits_{i=1}^nsumlimits_{j=1}^m(nmod i) imes(mmod j)=sumlimits_{i=1}^n(nmod i)sumlimits_{j=1}^m(mmod j)-sumlimits_{i=1}^n(nmod i)(mmod j))

设

(f(n,x)=sumlimits_{i=1}^ni(x/i))

则

(sumlimits_{i=1}^n(nmod i)=sumlimits_{i=1}^nn-i(n/i)=n^2-sumlimits_{i=1}^ni(n/i))

so

(egin{aligned}原式&=(n^2-f(n,n))(m^2-f(m,m))-sumlimits_{i=1}^n(nmod i)(mmod i)\&=(n^2-f(n,n))(m^2-f(m,m))-sumlimits_{i=1}^n(n-i(n/i))(m-i(m/i))end{aligned})

推导容斥的部分

(egin{aligned}sumlimits_{i=1}^n(n-i(n/i))(m-i(m/i))&=sumlimits_{i=1}^n(n-i(n/i))m-(n-i(n/i))i(m/i)\&=sumlimits_{i=1}^nnm-mi(n/i)-ni(n/i)+i^2(n/i)(m/i)\&=n^2m-msumlimits_{i=1}^ni(n/i)-nsumlimits_{i=1}^ni(m/i)+sumlimits_{i=1}^ni^2(n/i)(m/i)\&=n^2m-mf(n,n)-nf(n,m)+sumlimits_{i=1}^ni^2(n/i)(m/i)end{aligned})

最后一个 (sumlimits) 还是可以除数分块 (O(sqrt n)) 求。

so 预处理 (f(n,n),f(n,m),f(m,m)) 和上面的东西，就做完了。

最终复杂度 (O(sqrt n))

(Code:)

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
typedef long long ll;
const ll P=19940417;
inline ll read(){
	register ll x=0,f=0,ch=getchar();
	while('0'>ch||ch>'9')f^=ch=='-',ch=getchar();
	while('0'<=ch&&ch<='9')x=(x<<1)+(x<<3)+(ch^'0'),ch=getchar();
	return f?-x:x;
}
inline ll sum0(ll l,ll r){
	return 1LL*(r-l+1)*(l+r)%P*9970209%P;
}
inline ll sum(ll x){
	x%=P;
	return (1LL*x*(x+1)%P*((x+x)%P+1)%P*3323403%P)%P;
}
inline ll solve(ll t,ll x){
	register ll res=0;
	for(register int l=1,r;l<=t;l=r+1){
		if(x/l==0)r=t;
		else r=min(x/(x/l),t);
		res+=1LL*sum0(l,r)*(x/l)%P;
		res%=P;
	}	
	return res%P;
}
ll a,b,c;
ll n,m,f1,f2,f3,f4,ans;
inline void solved(){
	for(register int l=1,r;l<=n;l=r+1){
		r=min(n/(n/l),m/(m/l));
		f4+=1LL*(sum(r)%P-sum(l-1)%P+P)%P*(n/l)%P*(m/l)%P;
		f4%=P;
	}
}
signed main(){
	n=read(),m=read();
	if(n>m)n^=m^=n^=m;
	f1=solve(n,n);
	f2=solve(n,m);
	f3=solve(m,m);
	solved();
	ans=1LL*(1LL*n*n%P-f1+P)%P*(1LL*m*m%P-f3+P)%P;
	ans-=1LL*(1LL*n*n%P*m%P-(1LL*m*f1%P)-(1LL*n*f2%P)+f4)%P;
	ans=(ans+P)%P;
	printf("%lld
",ans);
	return 0;
}