2019ICPC西安邀请赛 - B. Product

打印的时候麻烦把：https://blog.csdn.net/skywalkert/article/details/50500009这个打印下来。

求$\prod\limits_{i=1}^{n} \prod\limits_{j=1}^{n} \prod\limits_{k=1}^{n} m^{gcd(i,j)[k|gcd(i,j)]} mod p$

欧拉定理：

$m^{ \sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} \sum\limits_{k=1}^{n} {gcd(i,j)[k|gcd(i,j)]} mod (p-1) }mod p$

算出指数直接套快速幂就可以了。考虑指数怎么算。为了方便把取模省略。

$\sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} \sum\limits_{k=1}^{n} {gcd(i,j)[k|gcd(i,j)]} $

这种带个整除的，一般都是交换求和到前面，然后后面算他的倍数。

$=\sum\limits_{k=1}^{n} \sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} {gcd(i,j)[k|gcd(i,j)]} $

带方括号的，枚举g，去括号。把g放在前面。

$=\sum\limits_{g=1}^{{n}g\sum\limits_{k=1}}{n} \sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} {[gcd(i,j)==g][k|g]} $

k与i,j无关了，提到前面。

$=\sum\limits_{g=1}^{{n}g\sum\limits_{k=1}}{n} [k|g] \sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} {[gcd(i,j)==g]} $

那么这个k就是d(g)了。

$=\sum\limits_{g=1}^{n}g*d(g) \sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} {[gcd(i,j)==g]} $

后面那个就把g提出来，好像还是很套路的。在什么鬼GCD统计里见过下式：

$\sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} {[gcd(i,j)g]}  = \sum\limits_{i=1}^{\lfloor\frac{n}{g}\rfloor} \sum\limits_{j=1}^{\lfloor\frac{n}{g}\rfloor} {[gcd(i,j)1]} $

那么这个两个数的互质对，规定i比j大于等于，那么就是欧拉函数，特别的，$\varphi(1)=1$。(i,j)和(j,i)互换重复计算(1,1)。

$ = \sum\limits_{i=1}^{\lfloor\frac{n}{g}\rfloor} \sum\limits_{j=1}^{\lfloor\frac{n}{g}\rfloor} {[gcd(i,j)==1]}  = 2 ( \sum\limits_{i=1}^{\lfloor\frac{n}{g}\rfloor} \varphi(i) ) - 1 $

这个东西用杜教筛板子都不用思考任何东西，直接弄出来（说起来我的杜教筛好像多个log）。

记上面的算式为 $S(\lfloor\frac{n}{g}\rfloor)$ ，那么原来的简化为：

$=\sum\limits_{g=1}^{n}g*d(g) S(\lfloor\frac{n}{g}\rfloor) $

换个喜欢的字母！！！

$=\sum\limits_{i=1}^{n}i*d(i) S(\lfloor\frac{n}{i}\rfloor) $

这些都跟g有关，太自闭了，大半时间卡在这里。后面的函数明显可以分块到 $O(\sqrt{n})$，对于一段连续的i求出结果。

要是可以计算前面的 $\sum\limits_{i=1}^{n}i*d(i)$ 的前缀和，那么可以再花$O(1)$求出一段连续的i的求和结果和后面相乘。

最后的问题就是怎么快速计算 $\sum\limits_{i=1}^{n}i*d(i)$ 这个狗东西。

首先可以线性预处理到5e6左右，花费一大堆内存，因为d(i)我是有线性筛的（虽然不会他的杜教筛）。

（菜鸡的）第一次重要的突破！

每个数的贡献和他的因子个数成线性。

然后我突发奇想，枚举每个因子d，枚举到 $\sqrt{n}$ ，每个因子会使得他的每个倍数都恰好多贡献1倍，那么就是 $d+2d+3d+...+\lfloor\frac{n}{d}\rfloor * d$。

等差数列求和嘛！

那么就变成：

$\sum\limits_{i=1}^{n}i*d(i)$

$=\sum\limits_{d=1}^{\lfloor\sqrt{n}\rfloor} \frac{d*(1+\lfloor\frac{n}{d}\rfloor)*(\lfloor\frac{n}{d}\rfloor)}{2}$

已经变成 $O(\sqrt{n})$求的了，当时出不了第三个样例，怎么想也是哦，能出就奇怪了。

（菜鸡的）第二次重要的突破！

233，看了半天才发现，上面不是也可以分块求的吗？对$\lfloor\frac{n}{d}\rfloor$分块啊。对一段连续的d，$\lfloor\frac{n}{d}\rfloor$都是同一个值，想了半天才领会！

那么就可以 $O(n^{\frac{1}{4}})$ 求出$\sum\limits_{i=1}^{n}i*d(i)$ 了，把这个狗东西记为 $T(i)$

$=\sum\limits_{i=1}^{n}T(i) S(\lfloor\frac{n}{i}\rfloor) $

杜教筛预处理 $O(n^{\frac{2}{3}})$ ，单次询问S(x)差不多 $O(1)$ ，虽然我的杜教筛比dq的多了一点东西。

单次询问T(i)差不多 $O(n^{\frac{1}{4}})$ 。意思是对单个i，处理出答案需要 $O(n^{\frac{1}{4}})$ 。

求和时，第一次分块对i分块用了 $O(\sqrt{n})$，每段连续的i共用一个T(i) S(\lfloor\frac{n}{i}\rfloor) ，所以就是 $O(n^{\frac{3}{4}})$ 。n是1e9，刚刚好够。

其他注意点：

1.应该用欧拉定理，这里卡了很久！幂取模不是直接取的，不要搞假算法！

2.涉及除法的运算，恰好求和公式可以约掉那个2，不要对那个数字取模！

总结：

$\sum\limits_{i=1}^{n}i*d(i)=\sum\limits_{d=1}^{\lfloor\sqrt{n}\rfloor} \frac{d*(1+\lfloor\frac{n}{d}\rfloor)*(\lfloor\frac{n}{d}\rfloor)}{2}$ 可以 $O(n^{\frac{1}{4}})$ 计算。

同理可得 $\sum\limits_{i=1}^{n}d(i)=\sum\limits_{d=1}^{\lfloor\sqrt{n}\rfloor} \frac{(1+\lfloor\frac{n}{d}\rfloor)*(\lfloor\frac{n}{d}\rfloor)}{2}$ 也是一样的。

说白了是我不会计算$d(i)$前缀和的方法，过于依赖线性筛了，去关注一下单个xx函数的求法吧。

当时写的代码，270分钟1A，好像挤到金银分隔线。

#include<bits/stdc++.h>

using namespace std;

typedef long long ll;

const int MAXN=5e6;

ll n,m,p,mod;

ll inv2;

ll qpow(ll x,ll n)

{

    ll res=1%p;

    while(n)

    {

        if(n&1)

            res=(res*x)%p;

        x=(x*x)%p;

        n>>=1;

    }

    return res%p;

}

ll Qarray[MAXN+5];

//Q(n)=\sum_i=1^n i*d(i)

inline ll Q(ll n)

{

    if(n<=MAXN)

        return Qarray[n];

    ll res=0;

    for(ll l=1,r;l<=n; l=r+1){

        r=n/(n/l);

        ll t=n/l;

        ll tmp1=((t+1)*t/2)%mod;

        ll tmp2=((l+r)*(r-l+1)/2)%mod;

        res=(res+(tmp1*tmp2%mod))%mod;

    }

    return res%mod;

}

unordered_map<int,ll> Sphi;

ll sump[MAXN+5];

int pri[MAXN+5],pritop;

bool notpri[MAXN+5];

ll d[MAXN+5],a[MAXN+5];

void sieve(int n=MAXN)

{

    pritop=0;

    notpri[1]=1;

    sump[1]=1;

    d[1]=1;

    for(int i=2; i<=n; i++)

    {

        if(!notpri[i])

        {

            pri[++pritop]=i;

            d[i]=2;

            a[i]=1;

            sump[i]=(i-1)%mod;

        }

        for(int j=1; j<=pritop&&i*pri[j]<=n; j++)

        {

            notpri[i*pri[j]]=1;

            if(i%pri[j])

            {

                d[i*pri[j]]=d[i]*d[pri[j]];

                a[i*pri[j]]=1;

                sump[i*pri[j]]=sump[i]*sump[pri[j]]%mod;

            }

            else

            {

                d[i*pri[j]]=d[i]/(a[i]+1)*(a[i]+2);

                a[i*pri[j]]=a[i]+1;

                sump[i*pri[j]]=sump[i]*pri[j]%mod;

                break;

            }

        }

    }

    Qarray[0]=0;

    Qarray[1]=1;

    for(int i=2; i<=n; i++)

    {

        sump[i]=(sump[i-1]+sump[i])%mod;

        Qarray[i]=(Qarray[i-1]+1ll*i*d[i])%mod;

    }

    return;

}

inline ll GetSphi(int n)

{

    if(n<=MAXN)

        return sump[n];

    if(Sphi.count(n))

        return Sphi[n];

    ll ret=n;

    if(ret%2==0)

        ret=(ret/2)*(1ll+n);

    else

        ret*=(1ll+n)/2;

    for(ll l=2,r;l<=n; l=r+1)

    {

        ll t=n/l;

        r=n/(t);

        ret-=(r-l+1)*GetSphi(t);

    }

    ret%=mod;

    return Sphi[n]=ret;

}

int main()

{

#ifdef local

    freopen("Yinku.in","r",stdin);

#endif // local

    scanf("%lld%lld%lld",&n,&m,&p);

    mod=p-1;

    sieve();

    ll ans=0;

    for(ll l=1,r; l<=n; l=r+1)

    {

        ll t=n/l;

        r=n/t;

        ll tmp=(Q(r)-Q(l-1)+mod)%mod;

        ll tmp2=GetSphi(t)%mod;

        ans=(ans+tmp*tmp2%mod)%mod;

    }

    //printf("ans=%lld\n",ans);

    ans=(ans*2ll)%mod;

    ans=(ans-Q(n)+mod)%mod;

    //printf("ans=%lld\n",ans);

    ll Ans=qpow(m,ans%mod);

    printf("%lld\n",Ans);

    return 0;

}

巴特西

2019ICPC西安邀请赛 - B. Product - 数论

其他注意点：

最新文章

热门文章