风暴止息

【题意简述】

有 $n$ 枚状态为 $0$ 或 $1$ 的定风楔。启动控制阵式 $i$ 会翻转编号为 $i$ 的全部正约数。每一步等概率随机启动一个阵式；若当前局面可以在不超过 $k$ 次操作内恢复，就立即改用最少操作方案。求总操作次数的期望乘以 $n!$ 后的值。

【Hint】

【提示】

从大到小可以唯一确定最少操作集合。随机操作只会使这个集合的大小增加或减少一，因此期望状态只需要记录集合大小。

【数据点 1∼51\sim 51∼5】

当 $n\leq 20$ 时，可以把全部状态压成二进制数，每个控制阵式对应一个异或掩码。对接管范围外的状态列期望方程，接管范围内直接使用最少操作数作为边界。

这个做法说明，真正需要研究的是当前局面距离全零状态还差多少次最优操作。

【数据点 6∼106\sim 106∼10】

从编号 $n$ 向 $1$ 扫描。若第 $i$ 枚定风楔失稳，只有阵式 $i$ 能改变它；之后编号更小的操作也不会再次影响位置 $i$。因此必须启动阵式 $i$，并翻转 $i$ 的全部正约数。

每一步选择都是被当前最大失稳位置强制确定的，所以得到的最少操作集合唯一。直接枚举约数可以在 $\mathcal{O}(n\sqrt n)$ 时间求出集合大小。

【数据点 11∼1511\sim 1511∼15】

设唯一最少操作集合大小为 $s$。随机选中集合内的阵式时，集合大小从 $s$ 变为 $s-1$；选中集合外阵式时，大小变为 $s+1$。因此向下转移的概率为 $\dfrac{s}{n}$，向上转移的概率为 $\dfrac{n-s}{n}$。

具体阵式编号已经不再影响随机过程，只需对状态数列期望方程。

【数据点 16∼1916\sim 1916∼19】

特殊性质保证初始最少操作数 $s\leq k$，随机阶段不会开始。答案就是

【正解】

令 $e_i$ 表示最少操作数从 $i$ 第一次变成 $i-1$ 所需随机操作次数的期望。显然 $e_n=1$。

当 $i<n$ 时，以 $\dfrac{i}{n}$ 的概率一步到达 $i-1$；以 $\dfrac{n-i}{n}$ 的概率先到达 $i+1$，随后需要 $e_{i+1}$ 回到 $i$，再重新等待 $e_i$。因此

把自环项移到左侧，得到

若初始最少操作数 $s\leq k$，答案就是 $s$。若 $s>k$，随机阶段依次从 $s$ 降到 $s-1,\ldots,k$，期望为

接管后还要执行 $k$ 次最优操作，所以总期望为

除法在模 $998{,}244{,}353$ 意义下乘逆元实现，最后再乘 $n!$。

最少操作集合可以从大到小枚举位置。正式 Std 在必须选择阵式 $i$ 时，枚举 $i$ 的全部正约数并翻转对应位置。单次枚举为 $\mathcal{O}(\sqrt i)$，最坏时间复杂度为 $\mathcal{O}(n\sqrt n)$。

也可以改为预处理每个编号的约数表，或从每个约数出发枚举倍数。后者的总枚举次数为调和级数

【复杂度分析】

正式 Std 求最少操作集合的最坏时间复杂度为 $\mathcal{O}(n\sqrt n)$，逐项快速幂求逆元为 $\mathcal{O}(n\log Mod)$。总时间复杂度为 $\mathcal{O}(n\sqrt n+n\log Mod)$，空间复杂度为 $\mathcal{O}(n)$。若使用调和级数维护与线性递推逆元，可以把时间复杂度优化为 $\mathcal{O}(n\log n)$。

【参考代码】

1/*2Author:Ehundategh3Date:2026/7/304Name:calm.cpp5You steal,I kill.6*/7#include <cmath>8#include <cstdio>9#include <cstring>10#include <algorithm>11#define MAXN 10001012using namespace std;13const int Mod=998244353;14inline int Add(int a,int b){return a+b>=Mod?a+b-Mod:a+b;}15inline int Del(int a,int b){return a-b<0?a-b+Mod:a-b;}16inline int Mul(int a,int b){return 1ll*a*b%Mod;}17int c,T,n,k,cnt,Fact,Dp[MAXN],Ans;18int Cond[MAXN];19inline int Inv(int a){20    int b=Mod-2,Ret=1;21    while(b){22        if(b&1) Ret=Mul(Ret,a);23        a=Mul(a,a);24        b>>=1;25    }26    return Ret;27}28 29void Solve(){30    scanf("%d%d",&n,&k);31    cnt=Ans=0;32    Fact=1;33    for(int i=1;i<=n;i++) scanf("%d",&Cond[i]),Fact=Mul(Fact,i);34    for(int i=n;i>=1;i--){35        if(Cond[i]){36            for(int j=1;j*j<=i;j++){37                if(i%j==0){38                    Cond[j]^=1;39                    if(j*j!=i) Cond[i/j]^=1;40                }41            }42            cnt++;43            Cond[i]=0;44        }45    }46    Dp[n]=1;47    if(cnt<=k){48        printf("%d\n",Mul(cnt,Fact));49        return;50    }51    for(int i=n-1;i>=1;i--) Dp[i]=Mul(Inv(i),Add(n,Mul(n-i,Dp[i+1])));52    for(int i=cnt;i>=k+1;i--) Ans=Add(Ans,Dp[i]);53    Ans=Mul(Add(Ans,k),Fact);54    printf("%d\n",Ans);55}56 57int main(){58    scanf("%d%d",&c,&T);59    while(T-->0) Solve();60    return 0;61}