洛谷P3768 简单的数学题

简述

　　这是道比较有意思的题。
　　首先用莫比乌斯反演推出答案：
　　令$g(x)=(\frac{x(x+1)}{2})^2$
　　

$$ans=\sum_{x=1}^nx^2g(\lfloor\frac{n}{x}\rfloor)\sum_{d|x}d\mu(\frac{x}{d})$$
　　后面那个卷积直接线性筛，然后直接做是 $O(n)$的，能过 $60$分。良心出题人部分分给的好足啊。
　　注意到一个事情，前面的那个 $g(\lfloor\frac{n}{x}\rfloor)$的取值是根号的，因此继续化下式子。
　　 $$ans=\sum_{x=1}^ng(\lfloor\frac{n}{x}\rfloor)x^2\sum_{d|x}d\mu(\frac{x}{d})$$
　　令 $$f(x)=x^2\sum_{d|x}d\mu(\frac{x}{d})$$
　　显然我们要预处理这个东西的前缀和才能做这道题，而且肯定是需要低于线性的算法（那就很可能是杜教筛了）。
　　这个东西很难看出什么门道，我们显然需要用狄利克雷卷积的运算律给他化简下。
　　令 $h(i)=i$
　　 $$f(x)=x^2(h*\mu)(x)$$
　　为了好看先把 $x^2$拿掉，我们看到里面有个 $\mu$，于是非常好奇地卷个 $1$
　　 $$h*\mu*1=h*(\mu*1)=h*\epsilon=h$$
　　咦？
　　也就是说 $(h*\mu)*1=h$
　　这不就是 $\varphi$吗， $\varphi*1=h$啊
　　太神奇了 （是我知道的太少了）
　　所以 $$f(x)=x^2\varphi(x)$$
　　这个东西怎么杜教筛啊，卷啥呀
　　冥思苦想了半天，还是搞出来了，令 $h_2(x)=x^2$
　　那么 $$(f*h_2)(x)=\sum_{d|x}f(x)h_2(\frac{x}{d})\\=\sum_{d|x}d^2\varphi(d)\frac{x^2}{d^2}\\=x^2\sum_{d|x}\varphi(x)=x^3$$

    　这样一个优美的结论着实令我惊讶。
    　　　　　　　　　　　　　　————引用自pty

　　令$s(n)=\sum_{i=1}^nf(x)$
　　则

$$h_2(1)s(n)=\sum_{i=1}^n(f*h_2)(i)-\sum_{i=2}^nh_2(i)s(\lfloor\frac{n}{i}\rfloor)\\s(n)=\sum_{i=1}^ni^3-\sum_{i=2}^nh_2(i)s(\lfloor\frac{n}{i}\rfloor)$$
　　然后你就不能挡我用杜教筛了。
　　你问我复杂度？
　　大概 $$O(N^{\frac{1}{2}}\times N^{\frac{1}{2}\times\frac{2}{3}})=O(N^\frac{5}{6})$$
　　确实比较玄学。

代码

#include <cstdio>
#include <algorithm>
#define maxn 4700000
#define ll long long
using namespace std;
ll N, p, prime[maxn], phi[maxn], f[maxn], s2[maxn], inv;
bool mark[maxn];
inline ll pow(ll a, ll b, ll p)
{
    ll ans=1, t=a;
    for(;b;b>>=1,t=t*t%p)if(b&1)ans=ans*t%p;
    return ans;
}
inline ll sqr(ll x){x%=p;return x*x%p;}
inline ll g(ll x){x%=p;return sqr((x*(1+x)>>1)%p);}
inline ll gets2(ll x)
{
    if(x<maxn)return s2[x];
    x%=p;
    return x*(x+1)%p*(2*x+1)%p*inv%p;
}
void init()
{
    ll i, j, d, x;
    phi[1]=1;
    for(i=2;i<maxn;i++)
    {
        if(!mark[i])prime[++*prime]=i, phi[i]=i-1;
        for(j=1;i*prime[j]<maxn;j++)
        {
            mark[i*prime[j]]=1;
            if(i%prime[j]==0){phi[i*prime[j]]=phi[i]*prime[j]%p;break;}
            phi[i*prime[j]]=phi[i]*phi[prime[j]]%p;
        }
    }
    for(i=1;i<maxn;i++)phi[i]=(phi[i]*sqr(i)%p+phi[i-1])%p;
    for(i=1;i<maxn;i++)s2[i]=(s2[i-1]+i*i)%p;
}
ll getf(ll n){return n<maxn?phi[n]:f[N/n];}
void calc(ll n)
{
    if(!n)return;
    ll i, last, t=N/n;
    if(n<maxn or f[t])return;
    f[t]=g(n);
    for(i=2;i<=n;i=last+1)
    {
        last=n/(n/i);
        calc(n/i);
        f[t]=(f[t]-(gets2(last)-gets2(i-1))*getf(n/i)%p)%p;
    }
}
void solve()
{
    ll ans=0, i, last, t;
    for(i=1;i<=N;i=last+1)
    {
        last=N/(N/i);
        calc(last);
        calc(i-1);
        ans=(ans+g(N/i)%p*((getf(last)-getf(i-1))%p)%p)%p;
    }
    printf("%lld",(ans+p)%p);
}
int main()
{
    ll ans=0, x;
    scanf("%lld%lld",&p,&N);
    inv=pow(6,p-2,p);
    init();
    solve();
    return 0;
}