zoukankan      html  css  js  c++  java
  • SDOI2015 约数个数和

    题目链接:戳我

    trick1——如何求约数个数和,变形

    [d(ij)=sum_{u|i}sum_{v|j}[gcd(u,v)=1] ]

    原式

    [=sum_{i=1}^Nsum_{j=1}^Msum_{u|i}sum_{v|j}[gcd(u,v)=1] ]

    [=sum_{u=1}^Nsum_{v=1}^Mlfloorfrac{N}{u} floorlfloorfrac{M}{v} floor[gcd(u,v)=1] ]

    [=sum_{i=1}^Nsum_{j=1}^Mlfloorfrac{N}{i} floorlfloorfrac{M}{j} floor[gcd(i,j)=1] ]

    (f(x)=sum_{i=1}^Nsum_{j=1}^Mlfloorfrac{N}{i} floorlfloorfrac{M}{j} floor[gcd(i,j)=x])

    [F(x)=sum_{x|d}f(d) ]

    [F(x)=sum_{x|d}sum_{i=1}^Nsum_{j=1}^Mlfloorfrac{N}{i} floorlfloorfrac{M}{j} floor[gcd(i,j)=d] ]

    trick2——消掉gcd

    [F(x)=sum_{i=1}^Nsum_{j=1}^Mlfloorfrac{N}{i} floorlfloorfrac{M}{j} floor[x|gcd(i,j)] ]

    [=F(x)=sum_{i=1}^{lfloorfrac{N}{x} floor}sum_{j=1}^{lfloorfrac{M}{x} floor}lfloorfrac{N}{ix} floorlfloorfrac{M}{jx} floor ]

    我们需要的是(f(1))

    [f(1)=sum_{d=1}^{min(N,M)}mu(d)F(d) ]

    [f(1)=sum_{d=1}^{min(N,M)}mu(d)sum_{i=1}^{lfloorfrac{N}{d} floor}lfloorfrac{N}{id} floorsum_{j=1}^{lfloorfrac{M}{d} floor}lfloorfrac{M}{jd} floor ]

    trick3——预处理(sum_{i=1}^{n}lfloorfrac{n}{i} floor)

    (g[i]=sum_{x=1}^ilfloorfrac{i}{x} floor)

    那么有

    [f(1)=sum_{d=1}^{min(n,m)}mu(d) imes g[n/d] imes g[m/d] ]

    直接数论分块即可。

    关于公式的证明:
    我没有严格公式化证明——
    我们可以发现,其实拼出来的每个数乘上他们的gcd,就是原先的因子,也就是互质形态相乘再乘上gcd的平方。现在需要解决的是这个值原先有没有被计算过,当然这个平方一个分配一个就是我们处理的原数,自然没有被计算过。而两个gcd放在一起乘上其中一个数,又因为我们的i,j的因子都是从小到大枚举的,所以大于当前数的数自然是没有被枚举计算过。
    然后严格证明捞了一位大佬的证明如下——

    代码如下:

    #include<iostream>
    #include<cstring>
    #include<cstdio>
    #include<algorithm>
    #define MAXN 50010
    int T,n,m,cnt;
    int mu[MAXN],sum[MAXN],not_prime[MAXN],prime[MAXN];
    long long g[MAXN];
    inline int min(int x,int y){return x>y?y:x;}
    inline void get_mu()
    {
        //not_prime[1]=1;
        mu[1]=1;
        for(int i=2;i<=50000;i++)
        {
            if(!not_prime[i]){prime[++cnt]=i;mu[i]=-1;}
            for(int j=1;j<=cnt&&prime[j]*i<=50000;j++)
            {
                not_prime[prime[j]*i]=1;
                if(i%prime[j]==0) break;
                else mu[i*prime[j]]=-mu[i];
            }
        }
        for(int i=1;i<=50000;i++) sum[i]=sum[i-1]+mu[i];
        for(int i=1;i<=50000;i++)
        {
            long long cur_ans=0;
            for(int l=1,r;l<=i;l=r+1)
                r=(i/(i/l)),cur_ans+=1ll*(r-l+1)*(i/l);
            g[i]=cur_ans;
        }
    }
    int main()
    {
        #ifndef ONLINE_JUDGE
        freopen("ce.in","r",stdin);
        #endif
        scanf("%d",&T);
        get_mu();
        //for(int i=1;i<=10;i++) printf("%d
    ",mu[i]);
        //for(int i=1;i<=10;i++) printf("%d
    ",sum[i]);
        while(T--)
        {
            int n,m;
            scanf("%d%d",&n,&m);
            long long ans=0;
            for(int l=1,r;l<=min(n,m);l=r+1)
            {
                r=min(n/(n/l),m/(m/l));
                ans+=1ll*(sum[r]-sum[l-1])*g[n/l]*g[m/l];
            }
            printf("%lld
    ",ans);
        }
    }
    
  • 相关阅读:
    HDU 2095 find your present (2) (异或)
    UESTC 486 Good Morning (水题+坑!)
    UVa 111 History Grading (简单DP,LIS或LCS)
    UVa 11292 Dragon of Loowater (水题,排序)
    HDU 1503 Advanced Fruits (LCS+DP+递归)
    UVa 10881 Piotr's Ants (等价变换)
    UVa 11178 Morley's Theorem (几何问题)
    HDU 1285 确定比赛名次(拓扑排序)
    .net Core的例子
    TCP与UDP的区别
  • 原文地址:https://www.cnblogs.com/fengxunling/p/10359896.html
Copyright © 2011-2022 走看看