• 51 NOD 1238 最小公倍数之和 V3


    原题链接

     

    最近被51NOD的数论题各种刷……(NOI快到了我在干什么啊!

    然后发现这题在网上找不到题解……那么既然A了就来骗一波访问量吧……

    (然而并不怎么会用什么公式编辑器,写得丑也凑合着看吧……

    $$
    ANS=sum_{i=1}^{n}sum_{j=1}^{n} frac{i*j}{gcd(i,j)}
    $$
    $$
    =sum_{d=1}^{n} d*sum_{i=1}^{leftlfloorfrac{n}{d} ight floor} sum_{j=1}^{leftlfloorfrac{n}{d} ight floor} i*j*[gcd(i,j)==1]
    $$
    $$
    =sum_{d=1}^{n} d*sum_{d'=1}^{leftlfloorfrac{n}{d} ight floor} f(leftlfloorfrac{n}{d*d'} ight floor)*d'^2*μ(d')
    $$
    $$
    =sum_{d'=1}^{n} d'^2*μ(d')*sum_{d=1}^{leftlfloorfrac{n}{d'} ight floor} d*f(leftlfloorfrac{n}{d*d'} ight floor)
    $$

    这里$f(x)=(frac{x*(x+1)}{2})^2$,即1到x中数字之间两两乘积之和。

     

    可以看出不同的$leftlfloorfrac{n}{d'} ight floor$只有根号n个,对于某一个$leftlfloorfrac{n}{d'} ight floor$,不同的$leftlfloorfrac{n}{d*d'} ight floor$也只有根号n个,那么把它们压在一起处理就好了。

    但是要做到这一点需要快速求出$d'^2*μ(d')$的前缀和,这时就要用上杜教筛了。

    记$g(x)=x^2*μ(x)$,$S(x)=sum_{i=1}^{x} g(i)$

    那么

    $  sum_{i=1}^{n} sum_{x|i} g(x)*(frac{i}{x})^2$

    $=sum_{i=1}^{n}  i*sum_{x|i} μ(x)$

    $=sum_{i=1}^{n} [i==1] $

    $= 1$

    同时

    $  sum_{i=1}^{n} sum_{x|i} g(x)*(frac{i}{x})^2$

    $=sum_{i=1}^{n} sum_{x=1}^{[frac{n}{i}]} g(x)*i^2$

    $=sum_{i=1}^{n} i^2*S[frac{n}{i}]$

    可以得到$S(n)=1-sum_{i=2}^{n} i^2*S[frac{n}{i}]$,预处理+哈希维护一波就行了

    但是不知道是我常数太大,还是方法没别人优越,卡了很久常数才卡过去。

    #include<ctime>
    #include<cstdio>
    #include<algorithm>
    #define ll long long
    #define N 1000001
    #define K 500000004
    using namespace std;
    ll read_p,read_ca,read_f;
    inline ll read(){
        read_p=0;read_ca=getchar();read_f=1;
        while(read_ca<'0'||read_ca>'9') {if (read_ca=='-') read_f=-1;read_ca=getchar();}
        while(read_ca>='0'&&read_ca<='9') read_p=read_p*10+read_ca-48,read_ca=getchar();
        return read_p*read_f;
    }
    const int HA=1e6+7;
    const int MOD=1e9+7;
    int MMH=0,p[N],num=0,mu[N],_ff[N],l[HA],Num,f[N],s[N],U=0;
    ll n,P;
    bool bo[N];
    struct na{int z,ne;ll y;}b[HA];
    inline void M(int &x){while (x>=MOD) x-=MOD;while(x<0)x+=MOD;}
    inline ll MM(ll x){return x<MOD&&x>-MOD?x:x%MOD;}
    inline void in(int x,ll y,int z){b[++Num].y=y;b[Num].z=z;b[Num].ne=l[x];l[x]=Num;}
    inline int _f(ll x){
        if (x>MOD) x%=MOD;
        if (x<N) return _ff[x];
        return MM(1LL*x*(x+1)>>1);
    }
    inline ll sqr(ll x){return x*x%MOD;}
    inline int F(ll x){
        int MMH=0;ll P;
        for (register ll i=1,j;i<=x;i=j+1) j=x/(P=x/i),M(MMH+=MM(sqr(_f(P))*(_f(j)-_f(i-1))));
        return MMH;
    }
    inline int Q(ll x){if (x>MOD) x%=MOD;return MM(MM(x*(x+1))*(x+x+1))*166666668%MOD;}
    inline int mmh(ll x){
        if (x<N) return s[x];
        register ll i,j;ll P;int o;
        for (i=l[o=x%HA];i;i=b[i].ne) if (b[i].y==x) return b[i].z;
        int MMH=0;
        for (i=2;i<=x;i=j+1) j=x/(P=x/i),M(MMH=MM(1LL*(Q(j)-Q(i-1))*mmh(P)+MMH));
        M(MMH=1-MMH);in(o,x,MMH);return MMH;
    }
    int main(){
        register int i,j;mu[1]=1;
        for (i=2;i<N;i++){
            if (!bo[i]) p[++num]=i,mu[i]=-1;
            for (j=1;j<=num&&i*p[j]<N;j++){
                bo[i*p[j]]=1;
                if (i%p[j]) mu[i*p[j]]=-mu[i];else {mu[i*p[j]]=0;break;}
            }
        }
        for (i=1;i<N;i++) M(f[i]=1LL*i*i*mu[i]%MOD+MOD),M(s[i]=s[i-1]+f[i]),_ff[i]=MM(1LL*i*(i+1)>>1);
        n=read();
        for (ll i=1,j;i<=n;i=j+1) j=n/(P=n/i),M(MMH+=MM(1LL*(mmh(j)-mmh(i-1))*F(P)));
        printf("%d
    ",MMH);
    }
    View Code

     

  • 相关阅读:
    6、linux中同步、互斥、阻塞(原子操作、信号量、阻塞)
    lightOJ-1199 Partitioning Game (SG函数)
    HDU-1013 Digital Roots
    HDU-1004 Let the Balloon Rise (STL map)
    HDU-1020 Encoding (字符串)
    POJ-2524 Ubiquitous Religions (并查集)
    POJ-1988 Cube Stacking (带权并查集)
    POJ-2236 Wireless Network (并查集)
    HDU-1002 A + B Problem II (模拟大数相加)
    HDU-1829 A Bug's Life (种类并查集)
  • 原文地址:https://www.cnblogs.com/Enceladus/p/5655267.html
Copyright © 2020-2023  润新知