好题啊！可惜没有链接！

题目描述

题解

计算方法和正解是一样的，但是推导过程用不着海森堡矩阵。

考虑把对角线上的元素全部替换为多项式 $C+x$，那么此矩阵的行列式就会变为一个多项式，可以代入 $x=1-C$ 求出答案。
组合一下选取的 $x$ 可以得到：
det ⁡ ( A ) = ∑ S ⊂ { 1 , 2 , . . . , n } det ⁡ ( B S ) x n − ∣ S ∣ \det(A)=\sum_{S\subset\{1,2,...,n\}}\det(B_S)x^{n-|S|} det(A)=S⊂{1,2,...,n}∑​det(BS​)xn−∣S∣其中 $B_S$ 是把 $A$ 的对角线换成 $C$ 然后限制在 $S$ 这个子集上的矩阵。比如：

观察 $B_S$ 可以发现，如果 $S$ 中的元素可以两两整除（也就是排序后，前一个整除后一个），那么 $B_S$ 是一个下三角矩阵，行列式为 $C^{|S|}$；
否则如果不能两两整除，那么一定有 $\det (B_S)=0$。这个可以递归证明：
如果 $S$ 中存在两个数不为其它任何数的因子，那么两个数所在的行都是全 $C$，那么矩阵一定不满秩；
反之，$S$ 中最大的数一定为所有数的倍数，那么它所在的列只有最后一行是 $C$，于是可以往删去最后一行一列的子矩阵递归下去。容易发现最后一定会出现上面不满秩的情况。

我们令 $r=\frac{C}{1-C}$，那么问题就转化为求满足元素两两整除的集合 $S$ 的贡献和：
$$Ans=(1-C{)}^n\sum _{S元素两两整除}{r}^{|S|}$$
令 g x = ∑ max ⁡ { i ∣ i ∈ S } = x r ∣ S ∣ g_x=\sum_{\max\{i|i\in S\}=x}r^{|S|} gx​=∑max{i∣i∈S}=x​r∣S∣（特别地，$g_1$ 表示 S = { 1 } S=\{1\} S={1} 和空集两种情况），那么可以枚举去掉最大值的集合，得到转移式子
$$g_i=r\sum _{j|i}{g}_j$$以及初值 $g_1=r+1$。
我们的目的是求出 $g_x$ 的前缀和，而这个只是看起来像狄利克雷卷积的式子不能用来 Min_25 筛，需要另找他法。

令 $s(n)={\sum }_{i=1}^n{g}_i$，我们其实可以直接枚举 $S$ 中的最大值和次大值的倍数关系得到状态转移：
$$s(n)=1+r+\sum _{i=2}^{n}s(\lfloor \frac{n}{i}\rfloor )$$这样我们就得到了一个直接根号分治枚举的做法。

然而这个做法复杂度为 $O(n^{\frac{3}{4}})$，只能过80分，还要优化。这里我们可以沿用杜教筛中的思想，先预处理出 $g_i$ 的前 $n^{\frac{2}{3}}$ 个值求出前缀和，然后剩下的部分根号分治枚举，这样就可以把根号分治的复杂度降至 $O(n^{\frac{2}{3}})$。

现在问题在于，怎么 $O(n)$ 线性地求出 $g_i$。考虑 $x$ 的唯一分解，我们发现若所有质因子的次数构成的可重集相同，那么 $g_x$ 相等（因为这本质上是一个往倍数走的路径计数，所以所有质因子等价）。我们线性筛求出所有元素的可重集的哈希值，那么哈希值相等的就只用求一次。

神奇的是，可重集的种类数非常少，少于 $\sqrt{线性筛范围}$，所以我们挑同一哈希值中最小的数直接枚举因子暴力求。

最后总复杂度即为 $O(n^{\frac{2}{3}})$。

代码

#include<bits/stdc++.h>//JZM yyds!!
#define ll long long
#define lll __int128
#define uns unsigned
#define fi first
#define se second
#define IF (it->fi)
#define IS (it->se)
#define END putchar('\n')
#define lowbit(x) ((x)&-(x))
#define inline jzmyyds
using namespace std;
const int MAXN=2e6+5;
const ll INF=1e18;
ll read(){
	ll x=0;bool f=1;char s=getchar();
	while((s<'0'||s>'9')&&s>0){if(s=='-')f^=1;s=getchar();}
	while(s>='0'&&s<='9')x=(x<<1)+(x<<3)+(s^48),s=getchar();
	return f?x:-x;
}
int ptf[50],lpt;
void print(ll x,char c='\n'){
	if(x<0)putchar('-'),x=-x;
	ptf[lpt=1]=x%10;
	while(x>9)x/=10,ptf[++lpt]=x%10;
	while(lpt>0)putchar(ptf[lpt--]^48);
	if(c>0)putchar(c);
}
const ll MOD=998244353;
ll ksm(ll a,ll b,ll mo){
	ll res=1;
	for(;b;b>>=1,a=a*a%mo)if(b&1)res=res*a%mo;
	return res;
}

unordered_map<ll,ll>mp;
mt19937 Rand(114514);
ll n,B,C,r,f[MAXN];
ll g[MAXN*10],mi[MAXN*10],hs[MAXN*10],ad[233];
bool nop[MAXN*10];
int pr[MAXN],le,m;
void sieve(int n){
	nop[0]=nop[1]=1;
	for(int a=2;a<=n;a++){
		if(!nop[a])pr[++le]=a,mi[a]=1,hs[a]=ad[1];
		for(int i=1,u;i<=le&&(u=pr[i]*a)<=n;i++){
			nop[u]=1;
			if(a%pr[i]==0){
				mi[u]=mi[a]+1,hs[u]=hs[a]-ad[mi[a]]+ad[mi[u]];
				break;
			}mi[u]=1,hs[u]=hs[a]+ad[1];
		}
	}g[1]=(r+1)%MOD;
	for(int x=2;x<=n;x++){
		ll&re=mp[hs[x]];
		if(!re){
			re=g[1];
			for(int i=2,lm=sqrt(x);i<=lm;i++)if(x%i==0){
				re+=g[i];
				if(x/i^i)re+=g[x/i];
			}re=re%MOD*r%MOD;
		}g[x]=re;
	}
	for(int i=2;i<=n;i++)(g[i]+=g[i-1])%=MOD;
}

int id(ll x){return x<=B?x:B+n/x;}
int main()
{
	freopen("bigben.in","r",stdin);
	freopen("bigben.out","w",stdout);
	n=read(),C=read(),B=sqrt(n),m=min(n,MAXN*10-50ll);
	// double MUDA=clock();
	if(!C)return print(1),0;
	if(C==1)return print(n>2?0:1),0;
	for(int i=1;i<=191;i++)ad[i]=Rand();
	r=C*ksm(MOD+1-C,MOD-2,MOD)%MOD,f[id(1)]=(r+1)%MOD;
	sieve(m);
	for(ll d=2,x;d<=n;d=x+1){
		if((x=n/(n/d))<=m){f[id(x)]=g[x];continue;}
		ll&res=f[id(x=n/(n/d))]=(r+1)%MOD;
		for(ll i=2,j,la;i<=x;i=la+1)
			j=x/i,la=x/j,(res+=(la-i+1)%MOD*r%MOD*f[id(j)])%=MOD;
	}
	print(f[id(n)]*ksm(MOD+1-C,n,MOD)%MOD);
	// printf("%.0f\n",clock()-MUDA);
	return 0;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88