题解：P10886 【MX-S3-T2】「FeOI Round 1」Journey

· · 题解

看题目有点像根号分治，但是观察数据范围，得出算法时间复杂度应该为 O(n)。显然不是根号分治。

再观察 g 的生成方式，比较一般，所以做法跟 g 的生成方式无关。

考虑计算 g_i 会对答案产生多少次贡献。当 a<i 时，对答案的贡献为 (n-i+1)\lfloor \frac{i-1}{c} \rfloor g_{i}；当 a=i 时，对答案的贡献为 (n-i+1)g_i；当 a>i 时，对答案没有贡献。所以答案为 \sum_{c=1}^{n}\sum_{i=1}^{n}(n-i+1)(\lfloor \frac{i-1}{c} \rfloor +1)g_i。这样就有一个 O(n^2) 的做法，可以拿 24 分。

考虑怎么优化，令 t=\lfloor \frac{i-1}{c} \rfloor，枚然后举 c,t，用类似于整除分块的思想，用差分进行区间修改。由于要要求 ct\le n，枚举的代码是调和级数，所以这个解法是 O(n\ln n) 的。核心代码如下：

for(int c=1;c<=n;c++)for(int t=1;c*t<=n;t++){
    b[c*t+1]+=t;if(c*(t+1)+1<=n)b[c*(t+1)+1]-=t;
}

这样不好优化，但是我们可以观察到这段代码和下面这段代码是等价的：

for(int c=1;c<=n;c++)for(int t=1;c*t<=n;t++){
    b[c*t+1]++;
}

至于为什么是等价的，留作习题答案略，读者自证不难。\ 其实就是把 t 分成多次计算贡献。

然后上面的代码得出的 b_i 就是 i-1 因数个数，可以线性筛 O(n) 处理。

代码如下：

#include<bits/stdc++.h>
#define ll long long
#define il inline
#define reg register
#define fmod(x) ((x)>mod?(x)-mod:(x))
using namespace std;

il void fread(ll &X){ll A=1;X=0;char C=getchar();while(!isdigit(C)&&C!='-')C=getchar();if(C=='-')A=-1,C=getchar();while(isdigit(C))X=(X<<1)+(X<<3)+(C^48),C=getchar();X*=A;}
il void fout(ll X){if(X<0)putchar('-'),X=-X;if(!X){putchar('0'),putchar(' ');return;}char C[25]{};int Len=0;while(X)C[++Len]=X%10+'0',X/=10;for(;Len;--Len)putchar(C[Len]);putchar(' ');}

const ll mod = 1e9+7,N=2e7+5;
ll n, A, B, C, g[N],ans=0;
ll sg,sgi,tmp,b[N],cnt=0;
bitset<N>v;int c[N],d[N],pri[N];

int main(){
    fread(n),fread(A),fread(B),fread(C),fread(g[n]);
    for(reg int i=n-1;i>=1;i--)g[i]=fmod(A*g[i+1]%mod*g[i+1]%mod+B*g[i+1]%mod+C);
    d[1]=1;for(reg ll i=2;i<=n;i++){
        if(!v[i]) pri[++cnt]=i,d[i]=2,c[i]=1;
        for(reg ll j=1;j<=cnt&&i*pri[j]<=n;j++){
            v[i*pri[j]]=1;if(i%pri[j]==0){
                c[i*pri[j]]=c[i]+1;d[i*pri[j]]=d[i]/(c[i]+1)*(c[i]+2);break;
            }d[i*pri[j]]=d[i]*2;c[i*pri[j]]=1;
        }
    }
    for(reg int i=1;i<=n;i++){
        sg=fmod(sg+g[i]),sgi=fmod(sgi+g[i]*i%mod);
        b[i]=d[i-1];b[i]=fmod(b[i]+b[i-1]);
        tmp=(n-i+1)*g[i]%mod;
        ans=fmod(ans+tmp*b[i]%mod);
    }ans=fmod(ans+n*(fmod(n*sg%mod+sg)-sgi+mod)%mod);fout(ans);
    return 0;
}