51nod1220约数之和

Time Limit 3s Memory Limit 128KB
Description
σ(x)表示x的约束和
求∑ni=1∑nj=1σ(ij)
Input
一个数n(0<n≤109)
Output
一个数答，对109+7取模。
Sample Input

1000

Sample Output

576513341

---

解题思路

莫比乌斯反演，杜教筛的各种技能
定义[x]的取值为1(x为真)否则为0

Ans=∑i=1n∑j=1n∑d|ijd

d|ij能不能化为一个d′|j呢，想到d|ij等同于dgcd(i,d)|j

d|ij要满足∃xy=d[x|i][y|j]，x取gcd(i,d)能供给完[x|i]也就是说，这时gcd(y,i)=1当且仅当y|j时d|ij成立，即d|ij↔dgcd(i,d)|j

∴

Ans=∑i=1n∑j=1n∑d=1ijd[dgcd(i,d)|j]=∑i=1n∑d=1ind⌊ngcd(i,d)d⌋=∑d′=1n∑i=1n∑d=1in[gcd(i,d)=d′]d⌊nd′d⌋      (d′表示gcd(i,d))

gcd(i,d)=d′?
考虑化为gcd(,)=1很简单啊，用id′∗d′代替i，用dd′∗d′代替d

Ans=∑d′=1n∑id′=1⌊nd′⌋∑dd′=1id′n[gcd(id′,dd′)=1]dd′∗d′⌊d′d∗n⌋

用i代替id′，用d代替dd′

Ans=∑d′=1n∑i=1⌊nd′⌋∑d=1in[gcd(i,d)=1]dd′⌊nd⌋=∑d′=1nd′∑i=1⌊nd′⌋∑d=1n[gcd(i,d)=1]d⌊nd⌋=∑i=1n∑d′=1⌊ni⌋d′∑d=1n[gcd(i,d)=1]d⌊nd⌋=∑i=1n12⌊ni⌋(⌊ni⌋+1)∑d=1n[gcd(i,d)=1]d⌊nd⌋

gcd(i,d)=1怎么处理呢，想到∑d|nμ(d)=[n=1]那么

[gcd(i,d)=1]=∑k|gcd(i,d)μ(k)=∑k[k|i][k|d]μ(k)

∴

Ans=∑i=1n12⌊ni⌋(⌊ni⌋+1)∑d=1nd⌊nd⌋∑d′=1n[d′|i][d′|d]μ(d′)=∑d′=1nμ(d′)∑id′=1⌊nd′⌋id′∗d′⌊nid′∗d′⌋∑dd′=1⌊nd′⌋12⌊ndd′∗d′⌋(⌊ndd′∗d′⌋+1)=∑d′=1nd′μ(d′)∑i=1⌊nd′⌋i⌊nid′⌋∑d=1⌊nd′⌋12⌊ndd′⌋(⌊ndd′⌋+1)

设

f(x)=xμ(x)

g(x)=∑i=1xi⌊xi⌋

h(x)=∑i=1x12⌊xi⌋(⌊xi⌋+1)

其实h和g是等价的

∑i=1x12⌊xi⌋(⌊xi⌋+1)=∑i=1x∑j=1⌊xi⌋j(i的出现次数和j的出现次数是一样的)=∑i=1x∑j=1⌊xi⌋i=∑i=1x⌊xi⌋i

∴

Ans=∑d=1nf(d)⋅g(⌊nd⌋)h(⌊nd⌋)=∑d=1ng(⌊nd⌋)2f(d)

我们想，⌊nd⌋只有2n√种

当i<=n√时，i∈[1,n√]
当i>n√时,ni<n√
∴⌊ni⌋只有2n√种
对于i，有j=max{x|⌊nx⌋=⌊ni⌋}=⌊n⌊ni⌋⌋

所以，分块处理
对于i，有∀j∈[⌊n⌊ni−1⌋⌋+1,⌊n⌊ni⌋⌋]⌊nj⌋=i
只需维护s(x)=∑xi=1f(x)即可，g(⌊ni⌋)直接算
有O(n)=∑⌊ni⌋O(ni−−√)=O(n34)

对于S(x)
∑d|nμ(d)=[n=1]
那么

1=∑i=1ni∑d|iμ(d)=∑d=1n∑i=1⌊ni⌋idμ(d)=∑d=1n∑i=1⌊ni⌋diμ(i)=∑d=1nd∑i=1⌊ni⌋f(i)=∑d=1nds(⌊nd⌋)

∴s(n)=1−∑d=2nds(⌊nd⌋)(将d=1项移项)

同样，分块处理∵⌊n4⌋=⌊⌊n2⌋2⌋，所以直接搜会有很多重复的项，哈希判重，O(n)=∑⌊ni⌋O(ni−−√)=o(n34)
预处理前n23项，复杂度可降到O(n23)

总复杂度O(n34)

#include<cstring>#include<cstdio>#include<cctype>#define N 3001001#define lim 43007#define mo 1000000007#define ny 500000004using namespace std;typedef long long ll;ll n,s[N],sum[lim],hash[lim],mu[N],pr[N];bool bz[N];int get(ll x){    int y=(int)(x%lim);    while(hash[y] && hash[y]!=x)y++,y=y==lim?0:y;return y;}ll S(ll a){    if(a<N)return s[a];    int x=get(a);    if(hash[x])return sum[x];    hash[x]=a;ll ans=1;    for(ll i=2,j;i<=a;i=j+1){        j=a/(a/i);        ans=(ans-(j+i)%mo*((j-i+1)%mo)%mo*ny%mo*S(a/i)%mo+mo)%mo;    }sum[x]=ans;return ans;}ll g(ll a){    ll ans=0;    for(ll i=1,j;i<=a;i=j+1){        j=a/(a/i);        ans=(ans+a/i%mo*((i+j)%mo*((j-i+1)%mo)%mo*ny%mo))%mo;    }return ans;}int main(){    s[1]=mu[1]=1;    for(ll i=2;i<N;i++){        if(!bz[i])mu[pr[++pr[0]]=i]=-1;        s[i]=(s[i-1]+mu[i]*i%mo+mo)%mo;        for(int j=1;j<=pr[0] && i*pr[j]<N;j++){            bz[i*pr[j]]=1;mu[i*pr[j]]=-mu[i];            if(i%pr[j]==0){mu[i*pr[j]]=0;break;}        }    }    scanf("%lld",&n);    ll ans=0;    for(ll i=1,j,las=0,now,t;i<=n;i=j+1,las=now){        j=n/(n/i);now=S(j);t=g(n/i);        ans=(ans+(now-las+mo)%mo*t%mo*t%mo)%mo;    }    printf("%lld",ans);}