学术篇51nod 1238 最小公倍数之和

Posted Rabbit House~❤

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了学术篇51nod 1238 最小公倍数之和相关的知识,希望对你有一定的参考价值。

这是一道杜教筛的入(du)门(liu)题目...


题目大意


\[ \sum_{i=1}^n\sum_{j=1}^nlcm(i,j) \]

一看就是辣鸡反演一类的题目, 那就化式子呗..
\[ \sum_{i=1}^n\sum_{j=1}^nlcm(i,j) \=\sum_{i=1}^n\sum_{j=1}^n\frac{ij}{gcd(i,j)} \=\sum_{i=1}^n\sum_{j=1}^n\sum_{k=1}^n\frac{ij}k[gcd(i,j)=k] \=\sum_{k=1}^n\sum_{i=1}^{\left\lfloor\frac nk\right\rfloor}\sum_{j=1}^{\left\lfloor\frac nk\right\rfloor}ijk[gcd(i,j)=1] \=\sum_{k=1}^nk\sum_{i=1}^{\left\lfloor\frac nk\right\rfloor}i\sum_{j=1}^{\left\lfloor\frac nk\right\rfloor}j[gcd(i,j)=1] \=\sum_{k=1}^nk\sum_{i=1}^{\left\lfloor\frac nk\right\rfloor}i\sum_{j=1}^{\left\lfloor\frac nk\right\rfloor}j\sum_{d=1}^{\left\lfloor\frac nk\right\rfloor}[d|i][d|j]\mu(d) \=\sum_{k=1}^nk\sum_{d=1}^{\left\lfloor\frac nk\right\rfloor}\mu(d) \sum_{i=1}^{\left\lfloor\frac nk\right\rfloor}[d|i]i \sum_{j=1}^{\left\lfloor\frac nk\right\rfloor}[d|j]j \设i=dq,j=dq, \原式=\sum_{k=1}^nk\sum_{d=1}^{\left\lfloor\frac nk\right\rfloor}\mu(d) \sum_{p=1}^{\left\lfloor\frac nk\right\rfloor}dp \sum_{q=1}^{\left\lfloor\frac nk\right\rfloor}dq \=\sum_{k=1}^nk\sum_{d=1}^{\left\lfloor\frac nk\right\rfloor}\mu(d)*d^2\sum_{p=1}^{\left\lfloor\frac nk\right\rfloor}p \sum_{q=1}^{\left\lfloor\frac nk\right\rfloor}q \\because\sum_{i=1}^ni=\frac{n(n+1)}2, \\therefore原式=\sum_{k=1}^nk\sum_{d=1}^{\left\lfloor\frac nk\right\rfloor}\mu(d)*d^2(\frac{{\left\lfloor\frac nk\right\rfloor}({\left\lfloor\frac nk\right\rfloor}+1)}{2})^2 \设t=kd, \原式=\sum_{t=1}^n(\frac{{\left\lfloor\frac nt\right\rfloor}({\left\lfloor\frac nt\right\rfloor}+1)}{2})^2\sum_{d|t}\mu(d)*td \设g(x)=\sum_{d|x}\mu(x)*xd, \则原式=\sum_{t=1}^n(\frac{{\left\lfloor\frac nt\right\rfloor}({\left\lfloor\frac nt\right\rfloor}+1)}{2})^2g(t) \]
然后前面这一堆我们可以分块\(O(1)\)搞出来, 所以问题的关键就变成了\(g(x)\)的前缀和.
因为题目中\(n\leq10^{10}\), 所以这里要用\(O(n^\frac23)\)的杜教筛.
但是杜教筛的框架我们是知道的,
\[ 设要求前缀和的函数为f(x), 其前缀和为s_f(x) \设有好求前缀和的积性函数g(x)和h(x), 使得(f*g)(x)=h(x), \设h(x)的前缀和为s_h(x), 则有 \s_f(x)=\frac{s_h(x)-\sum_{d=2}^ns_f({\left\lfloor\frac nd\right\rfloor})g(d)}{g(1)} \]
所以我们的任务又变成了找到合适的\(g(x)\)\(h(x)\).
在某篇歪果仁的blog上, 他说:

With a little luck, we can notice that Id(l)?=?g?*?Id2(l),

  • 注: \(id(x)=x\)(以下称为\(n\))\(id2(x)=x^2\)(以下称为\(n_2\))

But 这运气也太好了点吧? 我怎么就找不到这种函数呢?!
不过我们还是要考虑一下的... (这波化式子的时候注意乘号\(\cdot\)和卷积号\(*\))
\[ g(x)=\sum_{d|x}\mu(d)xd \=\sum_{d|x}x^2\mu(d)\frac dx \=x^2\sum_{d|x}\mu(d)\frac dx \=x^2\cdot(\mu*n) \(g*n_2)(x) \=(x^2\cdot\mu*n)*n_2 \cdots ① \积性函数有一个性质, 如果两边的乘积都含有某一项的时候, 可以把这一项提出来. \所以我们可以把x^2都提出来, (证明可以用定义式推咯~) \①=(x^2\cdot\mu*n_2)*n \=x^2\cdot(\mu*1)*n \=x^2\cdot\epsilon*n =n \]
所以我们就证明出来了\((g*n_2)(x)=n(x)\)
这样就可以直接杜教筛了= =
然后就是注意细节了OvO
比如非常坑人的\(10^{10}*(1e9+7)\)会爆long long!!
所以要开unsigned...
就这样吧= =
代码:

#include <map>
#include <cstring>
#include <iostream>
using namespace std;
typedef unsigned long long LL;
#define ri register int
#define rll register unsigned long long 
const int p=1e9+7;
int qpow(int a,int b,int s=1){
    for(;b;b>>=1,a=1LL*a*a%p)
        if(b&1) s=1LL*s*a%p;
    return s;
}
const int inv2=qpow(2,p-2),inv4=qpow(4,p-2),inv6=qpow(6,p-2);
const int N=4800000,_N=N+10;
const int SQ=233333;
LL n;
map<LL,int> mp;
map<LL,int>::iterator it;
int g[_N],pr[N/2],tot;
bool notp[_N];
void shai(){
    g[1]=1; ri i,j; rll k;
    for(i=2;i<=N;++i){
        if(!notp[i])
            pr[++tot]=i,g[i]=(-1LL*i*(i-1)%p+p)%p;
        for(j=1;j<=tot&&(k=1LL*i*pr[j])<=N;++j){
            notp[k]=1;
            if(i%pr[j]==0){g[k]=1LL*g[i]*pr[j]%p;break;}
            g[k]=1LL*g[i]*g[pr[j]]%p;
        }
    }
    for(i=2;i<=N;++i) g[i]=(g[i]+g[i-1])%p;
}
inline int calcsqr(LL x){x%=p;return x*(x+1)%p*(2*x+1)%p*inv6%p;}
inline int calcsum(LL x){x%=p;return x*x%p*(x+1)%p*(x+1)%p*inv4%p;}
int calc(LL x){
    if(x<=N) return g[x];
    it=mp.find(x);
    if(it!=mp.end()) return it->second;
    int ans=1LL*x%p*(x+1)%p*inv2%p;
    LL last=0;
    for(rll i=2;i<=x;i=last+1){
        last=x/(x/i);
        ans=(ans-1LL*calc(x/i)*(calcsqr(last)-calcsqr(i-1)+p)%p+p)%p;
    }
    return mp[x]=(ans%p+p)%p;
}
int main(){
    shai();cin>>n;
    LL last=0; int ans=0;
    for(rll i=1;i<=n;i=last+1){
        last=n/(n/i);
        ans=(ans+1LL*calcsum(n/i)*(calc(last)-calc(i-1)+p)%p)%p;
    }
    cout<<ans;
}

这破题我连推带写加总结写了整整一天!!!
窝还是太弱了..

以上是关于学术篇51nod 1238 最小公倍数之和的主要内容,如果未能解决你的问题,请参考以下文章

51nod1238 最小公倍数之和 V3

[51nod1238] 最小公倍数之和 V3(杜教筛)

51nod1238 最小公倍数之和 V3

51nod1238 最小公倍数之和 V3

51nod - 1363 - 最小公倍数之和 - 数论

51nod-1363: 最小公倍数之和