为什么你现在叫 CyanSineSin 了？？

· · 题解

ZJOI2022 Day2 的题都是非人力可及牛逼题目吗？好吓人啊！！！

先随便讲点暴力，因为这显然是线性变换，加上一些聊胜于无的特判就可以拿 50。

注意到有个奇怪的特殊限制，一个是 n=2^k，一个是 n=98304 = 2^{15} \cdot 3。我们对这俩东西研究一下。

首先是 n=2^k。拉一次拉面，会发现 a_1=a_2,a_3=a_4 ,\cdots ,a_{n-1}=a_n；再拉一次，会发现 a_1=a_2=a_3=a_4 ,\cdots ,a_{n-3}=a_{n-2}=a_{n-1}=a_n……最终所有数会变成一样的，并不再改变。

再考虑 n=2^{15} \cdot 3。类似的，在第十五次操作时，会分成三个极长段，每段内数相同，长度均为 2^{15}。再操作一次就变成了前面一个 2^{16} 的极长段后面一个 2^{15} 的极长段。继续操作，权值仍然会改变，但是这样的极长段的样式不会改变。

经过多次手玩，记最大的 x 满足 n \equiv 0 \pmod{2^x} 为 t，在 t+1 次操作后会形成如下的形状：

a_1^{2^{t+1}}a_2^{2^{t+1}}a_3^{2^{t+1}} \cdots a_{m-1}^{2^{t+1}} a_m^{2^t}

最后一段长度只有前面的一半。注意到形成这个形状的最大操作次数是 O(\log n) 级别的，我们先预处理掉前 O(\log n) 次操作，然后考虑把现在的这个 a 序列抽出来，研究其变化。

那么操作里面有个非常多余的东西就是除以二，这个除以二什么时候做都可以，可以从线性变换的角度考虑。那我们在最后除以 2^k 即可（要特判 n=2^p）。

假设现在的 a 序列是 a_1,a_2,\cdots a_m 并记 mid = \lceil \frac{m}{2} \rceil，那么操作一次之后是 a_1+a_m,a_1+a_{m-1},a_2+a_{m-1},a_2+a_{m-2},\cdots ,a_{mid-1}+a_{mid},2a_{mid}。

这个形式显得有规律可循，考虑差分。记之前的差分数组为 b_1,b_2,\cdots ,b_{m-1}，那么现在的差分数组是 -b_{m-1},b_1,-b_{m-2},b_2 \cdots ,b_{mid-1}。这是这个题离谱的第一点。

现在可以做到 O(n^2+q) 了，还可以更优秀吗？

这个数组的变换是非常有规律的，但是出现了 -1 这个系数显得不优美，我们考虑往差分数组后面补一个 -b_{m-1},-b_{m-2},\cdots ,-b_{1}，然后这是一个置换的形式。现在每次询问还原序列就好了，时间复杂度 O(nq)。

然后补充一个离谱的观察。这个置换拆成若干个置换环，满足大小最大的置换环是其他置换环大小的倍数……因此可以保证其循环节在 O(n) 范围内。进一步的，最大的置换环大小等于 2 在 \bmod (n+1) 下乘法的阶（意即最大的置换环大小 x 满足 2^x \equiv 2 \pmod{(n+1)} 且 x 最小）。这是一个非常牛逼的观察……离谱的第二步。实则有迹可循，你可以发现在上面那个差分数组 d 中，满足 d_{2i \bmod (m+1)}' = d_i。这个观察比较重要。

还能够优秀！考虑通过差分数组快速计算 x 位置上的值。要做的两步是，一个差分值在算 a_x 的时候造成了多少次贡献，以及 a_1 的值是多少。

首先第一个观察是，因为我们将除以二这个东西提到了最后，每次操作数列的权值和是不变的。那么求 a_1 并不是一件困难的事情，其表示成总和减去差分数组带上一些系数，最后除以 n。然后算 x 的话，也是近似于差分数组带系数。不妨记这个系数为 p_i 与 q_i，其乘上新的差分数组 d_i^{(k)} 才有意义。同时这个 k 优化后是 O(n) 级别的。

记 S 为变换前 a 的和，那么有：

a_1 \times n + \sum_{i=1}^{m-1}d_i^{(k)}p_i = S a_x = a_1 + \sum_{i=1}^{m-1}d_i^{(k)}q_i

（容易发现 q 后面一段可以是 0，同时 p,q 是容易得出的。）

我们要对所有的 k 全部求出一个答案（同时 k 最大是 O(n) 级别的）。相当于求 \displaystyle \sum_{i=1}^{m-1}d_ip_i, \sum_{i=1}^{m-1} d_i^{(1)}p_i,\sum_{i=1}^{m-1} d_i^{(2)}p_i,\cdots（其中的 p 显然换成 q 也要算一次）。

那么我们将 d 表示成置换的形式。但是变换的形式就是 d_{2i \bmod (m+1)}' = d_i 就显得比较麻烦啊。那就把置换环全部拆出来，这样的话这个麻烦的变换形式就显得简单一些，类似于循环卷积。

至于循环卷积，在这个题中举一个例子出来。比如一个长为 P 的置换环 d，置换 k 次后对上一个系数数组 p，有：

\sum_{i=0}^{P-1} d_i p_{i+k \bmod P}

同时对每个 k 都要求一个这个。这里就显得非常的简单，把 p 翻转做循环卷积就好了。可以用 NTT 实现。

注意到不同的置换环，如果大小相同造成贡献的位置也相同，可以提到一起做。这样的话大小不同的置换环个数只有 O(\sqrt n) 级别，可以做到 O(n\log n+q \sqrt n)。

还有一点就是我们一直在讨论上面那个最后一段长度是前面的段的一半的 case，因此如果 n=2^k 需要特判。

当然根据我们上面那个比较离谱的观察，因为最大的置换环大小等于 2 在 \bmod (n+1) 下乘法的阶，且其他置换环大小是最大的置换环的因数，因此整个 a 会产生循环节，且长度非常有限。这样的话记最大的置换环长度为 P，直接预处理所有 k\bmod P=0,1,\cdots P-1 的答案就好了。这样就可以做到 O(n \log n + n d(n) + q) 了……吗？

注意到 k 上面的性质都是我们在 a 开始形成上面那种特殊形状后去掉除以二的操作得到的，我们还要求一个 \dfrac{1}{2^k}，这个 k 是没有经过循环节取模优化的真实的 k，范围仍然是 10^{18}……所以先欧拉定理降次，然后光速幂就好了。

细节比较多。完整代码戳这里。

int a[2000005],test;
int hst[25][2000005],b[2000005],d[2000005],m;
int cir[2000005];
bool vis[2000005];
int cbid[2000005];
int fcir[2000005];
void Solve()
{
    int n, q, x;
    long long k_max;
    scanf("%d%d%d%lld", &n, &q, &x, &k_max);
    int Sum=0;
    for (int i = 1; i <= n; i ++)   scanf("%d",&a[i]),Sum=Add(Sum,a[i]);
    for(int i=1;i<=n;++i)   hst[0][i]=a[i];
    for(int i=1;i<=24;++i)
    {
        for(int j=1;j<=n;++j)   b[j]=Add(hst[i-1][j],hst[i-1][n-j+1]);
        for(int j=1;j<=n;++j)   hst[i][j]=Mul(b[(j+1)/2],inv2);
    }
    if((n&(-n))==n)
    {
        LL ans=0;
        for(int C=1;C<=q;++C)
        {
            LL k=rd(seed)%k_max;
            k=min(k,24ll);
            ans^=(LL(C)*hst[k][x]);
        }
        printf("%lld\n",ans);
        return ;
    }
    int t=0;
    {
    int tmp=n;
    while(tmp%2==0)
    {
        tmp>>=1;
        ++t;
    }
    ++t;
    }
    int tp=t;
    {
    int pm=0;
    for(int i=1;i<=n;i+=(1<<t)) b[++pm]=hst[t][i];
    t=1<<t;
    m=0;
    for(int i=2;i<=pm;++i)  d[++m]=Sub(b[i],b[i-1]);
    for(int i=m;i;--i)  d[++m]=MOD-d[i];
    }
    for(int i=1;i<=m;++i)   cir[2*i%(m+1)]=i;
    for(int i=1;i<=m;++i)   vis[i]=false;
    vector<Poly> Cir;
    for(int i=1;i<=m;++i)
    {
        if(vis[i])  continue;
        int u=i;
        Cir.push_back({});
        while(!vis[u])
        {
            vis[u]=true;
            Cir.back().push_back(u);
            u=cir[u];
        }
    }
    int K=1;
    for(auto st:Cir)    K=max(len(st),K);
    for(int i=0;i<=K;++i)   cbid[i]=fcir[i]=0;
    int inv=QuickPow(n);
    vector<Poly> CirBuc;
    for(auto st:Cir)
    {
        Poly A,B;
        int len=len(st);
        A.resize(len),B.resize(len);
        for(int i=0;i<len;++i)  A[i]=d[st[i]];
        for(int i=0;i<len;++i)  B[i]=MOD-Mul(inv,st[i]<=m/2?Sub(n,Mul(st[i],t)):0);
        for(int i=0;i<len;++i)  B[i]=Add(B[i],int(st[i]*t+1<=x));
        reverse(B.begin(),B.end());
        Poly C=A*B;
        for(int i=0;i<len-1;++i)    C[i]=Add(C[i],C[i+len]);
        if(!cbid[len])
        {
            CirBuc.push_back({});
            CirBuc.back().resize(len);
            cbid[len]=len(CirBuc);
        }
        int id=cbid[len]-1;
        for(int i=0;i<len;++i)  CirBuc[id][i]=Add(CirBuc[id][i],C[i]);
    }
    for(int len=1;len<=K;++len)
    {
        if(cbid[len])
        {
            int id=cbid[len]-1;
            for(int i=0,j=len-1;i<K;++i,j=(j+1)%len)    fcir[i]=Add(fcir[i],CirBuc[id][j]);
        }
    }
    LL ans=0;
    for(int C=1;C<=q;++C)
    {
        LL k=rd(seed)%k_max;
        int ret=0;
        if(k<=tp)   ret=hst[k][x];
        else
        {
            k-=tp;
            int coe=cp.Inv(k);
            k%=K;
            int s=Mul(fcir[k],coe);
            ret=Add(s,Mul(inv,Sum));
        }
        ans^=LL(C)*LL(ret);
    }
    printf("%lld\n",ans);
}