【51nod】1238 最小公倍数之和 V3 杜教筛

【题意】给定n,求Σi=1~nΣj=1~n lcm(i,j),n<=10^10。

【算法】杜教筛

【题解】就因为写了这个非常规写法,我折腾了3天……

$$ans=\sum_{i=1}^{n}\sum_{j=1}^{n}lcm(i,j)$$

$$g(n)=n*\sum_{i=1}^{n}\frac{i}{(n,i)}$$

那么

$$ans(n)=2*g(n)-\sum_{i=1}^{n}i$$

枚举gcd,化简g(n)。

$$g(n)=n*\sum_{d|n}1/d\sum_{i=1}^{n}i*[(n,i)=d]$$

令i=i/d

$$g(n)=n*\sum_{d|n}1/d\sum_{i=1}^{n/d}id*[(n/d,i)=1]$$

$$g(n)=n*\sum_{d|n}\sum_{i=1}^{n/d}i*[(n/d,i)=1]$$

由于

$$\sum_{i=1}^{n}[(n,i)=1]*i=\frac{n*\varphi(n)+[n==1]}{2}$$

所以代入,得

$$g(n)=n*\sum_{d|n}\frac{d*\varphi(d)+[d==1]}{2}$$

这里需要注意取整的问题,当d>1时d*φ(d)一定是偶数,当d=1时d*φ(d)=1就必须结合[d==1],于是可以化简成下面的形式。

$$g(n)=\frac{1}{2}n(1+\sum_{d|n}d*\varphi(d))$$

$$g(n)=\frac{1}{2}(n+n*\sum_{d|n}d*\varphi(d))$$

将上式代入ans,得

$$ans=2*\sum_{i=1}^{n}\frac{1}{2}(i+i*\sum_{d|i}d*\varphi(d)))-\sum_{i=1}^{n}i$$

$$ans=\sum_{i=1}^{n}i*\sum_{d|i}d*\varphi(d)$$

$$f(n)=n*\sum_{d|n}d*\varphi(d)$$

那么

$$ans=F(n)=\sum_{i=1}^{n}f(n)$$

 

★呼,经过上面一系列的化简,我们终于来到了杜教筛——求f(n)的前缀和。

 

$$f(n)=n*\sum_{d|n}d*\varphi(d)$$

将f表示成狄利克雷卷积的形式,根据点积的卷积分配律(乱起的名字)。

$$f=id \cdot (1*(id \cdot \varphi))=id*(id^2 \cdot \varphi)$$

发现其中有id^2的形式,我们知道同阶幂函数卷积有奇效www。(实际上应该是相同完全积性函数卷积有奇效)

$$g=id^2$$

$$f*g=(id^2 \cdot \varphi)*id*id^2=[(id^2 \cdot \varphi)*id^2]*id$$

双重卷积非常麻烦,考虑先化简方括号内的卷积

$$([(id^2 \cdot \varphi)*id^2])(i)=\sum_{d|n}d^2*\varphi(d)*\frac{n^2}{d^2}=n^2\sum_{d|n}\varphi(d)=n^3$$

成功化简!我们可以把f*g表示出来了!

$$(f*g)(i)=\sum_{d|n}d^3*\frac{n}{d}=n*\sum_{d|n}d^2$$

jiry_2在她的博客中表示这个柿子的前缀和是非常好求的。

我:???

假设h=f*g,那么

$$H(i)=\sum_{i=1}^{n}i*\sum_{d|i}d^2=\sum_{i=1}^{n}i^2\sum_{d=1}^{n/i}d*i=\sum_{i=1}^{n}i^3*\frac{\frac{n}{i}*(\frac{n}{i}+1)}{2}$$

然后就可以进行分块取值优化了,注意G和H的求解复杂度为O(1)~O(√n)对杜教筛的总复杂度都没有影响。

最终杜教筛的形式是

$$F(n)=H(n)-\sum_{i=2}^{n}i^2*F(\frac{n}{i})$$

 

#include<cstdio>
#include<cstring>
#include<algorithm>
#define int long long
using namespace std;
const int MOD=1e9+7,preN=10000000,v=(MOD+1)/2;
int N,f[preN+5],phi[preN+5],e[preN+5],a[200010],prime[1000010],tot;
int A(int n){return n%MOD*(n%MOD+1)%MOD*v%MOD;}
int B(int n){return n%MOD*(n%MOD+1)%MOD*(2*n%MOD+1)%MOD*((MOD+1)/6)%MOD;}
int K(int x){return x*x%MOD;}
int C(int n){return K(A(n));}
int query(int n){
    int pos=1,ans=0;
    for(int i=1;i<=n;i=pos+1){
        pos=n/(n/i);
        ans=(ans+(C(pos)-C(i-1)+MOD)*A(n/i)%MOD)%MOD;
    }
    return ans;
}
int solve(int n){
    if(n<=preN)return f[n];
    if(~a[N/n])return a[N/n];
    int ans=query(n);
    int pos=2;
    for(int i=2;i<=n;i=pos+1){
        pos=n/(n/i);
        ans=(ans-(B(pos)-B(i-1)+MOD)*solve(n/i)%MOD+MOD)%MOD;
    }
    return a[N/n]=ans;
}
#undef int
int main(){
#define int long long
    scanf("%lld",&N);
    phi[1]=1;f[1]=1;e[1]=1;
    for(int i=2;i<=preN;i++){
        if(!e[i]){
            phi[prime[++tot]=i]=i-1;
            e[i]=i;
            f[i]=(i*phi[i]+1)%MOD;
        }
        for(int j=1;j<=tot&&i*prime[j]<=preN;j++){
            int k=i*prime[j];
            if(i%prime[j]==0){
                e[k]=e[i]*prime[j];
                phi[k]=phi[i]*prime[j];
                f[k]=(f[i]+f[i/e[i]]*phi[e[k]]%MOD*e[k]%MOD)%MOD;
                break;
            }
            e[k]=prime[j];
            phi[k]=phi[i]*(prime[j]-1);
            f[k]=f[i]*f[prime[j]]%MOD;
        }
    }
    for(int i=2;i<=preN;i++)f[i]=(f[i-1]+f[i]*i%MOD)%MOD;
    memset(a,-1,sizeof(a));
    printf("%lld",solve(N));
    return 0;
}
View Code

 

 

 

 

1.F(n)的预处理:在i%prime[j]时,f[i*prime[j]]=f[i]+f[i/e[i]]*φ(e[i]*prime[j])*(e[i]*prime[j]),其中e[i]是 i 的所有最小素因子乘积(p1^k1)。

2. 1~3次幂前缀和

$$\sum_{i=1}^{n}i=\frac{n(n+1)}{2}$$

$$\sum_{i=1}^{n}i^2=\frac{n(n+1)(2n+1)}{6}$$

$$\sum_{i=1}^{n}i^3=(\frac{n(n+1)}{2})^2$$

3.x的逆元是(MOD+1)/x,如果这是个整数。

4.全程long long的话,define比较好。而且输入的n是long long的话,n*n会爆long long。

 

 

 

【另一种写法】

 

$$ans=\sum_{i=1}^{n}\sum_{i=1}^{n}\frac{i*j}{(i,j)}$$

 

直接枚举gcd值

 

$$ans=\sum_{d=1}^{n}1/d\sum_{i=1}^{n/d}\sum_{i=1}^{n/d}d^2*i*j*[(i,j)=1]$$

 

$$ans=\sum_{d=1}^{n}d\sum_{i=1}^{n/d}\sum_{j=1}^{n/d}i*j*[(i,j)=1]$$

 

这里看到[(i,j)=1]很容易想到莫比乌斯反演,但是这个问题如果反演会变得相当复杂,应该留到后面用n*φ(n)/2化简。

 

从上式已经可以看出分块取值优化的形式了。

 

 

$$s(n)=\sum_{i=1}^{n}\sum_{j=1}^{n}i*j*[(i,j)=1]$$

 

那么

 

$$ans=\sum_{d=1}^{n}d*s(n/d)$$

 

为了将s(n)表示成前缀和的形式,将矩形转化为上三角。

 

$$s(n)=2*(\sum_{i=1}^{n}\sum_{j=1}^{i}i*j*[(i,j)=1])-1$$

 

$$s(n)=2*P(n)-1$$

 

$$p(n)=n*\sum_{i=1}^{n}i*[(n,i)=1]$$

 

然后就可以转化了

 

$$p(n)=n*\frac{n*\varphi(n)+[n==1]}{2}=\frac{1}{2}(n^2\varphi(n)+n)$$

 

 

$$f(n)=n^2\varphi(n)$$

 

那么s(n)可以表示为

 

$$s(n)=F(n)+\frac{n(n+1)}{2}-1$$

 

然后我们就可以用杜教筛求解f(n)的前缀和。

 

$$f=id^2 \cdot \varphi$$

 

$$g=id^2$$

 

$$f*g=id^3$$

 

因为分块和杜教筛都是求n/i,所以复杂度并列,最终O(n^2/3)。

 

posted @ 2018-02-28 16:26  ONION_CYC  阅读(195)  评论(0编辑  收藏