【XSY4350】摆(行列式,数论,杜教筛)

题面

题解

首先我们将原矩阵写成 \(A+B\),其中 \(B\) 全是 \(C\),那么 \(A\) 的第 \(i\) 行就只有其倍数处有值,且 \(A_{i,i}=1-C,A_{i,j(i|j\land i\neq j)}=-C\)

那么原来的行列式就变成了:

\[\sum_{p}(-1)^{\operatorname{sgn}(p)}\sum_{S\subseteq \{1,\cdots,n\}}\prod_{i\in S}A_{i,p_i}\prod_{i\not\in S}B_{i,p_i} \]

考虑在 \(A\) 中选的点数(\(|S|\)),发现它不可能小于 \(n-1\),否则对于 \(B\) 要求的就是一个大小大于 \(1\)、全是 \(C\) 的矩阵的行列式,它恒为 \(0\)

那么假设 \(B\) 中选的那个点(\(B\) 中没选点的情况很好算,略)为 \((i,p_i)\)。考虑把行列式问题放在图论上(注意到 \((-1)^{\operatorname{sgn}(p)}\)\((-1)^{n-\text{环个数}}\),所以系数的问题也是可以处理的),相当于说 \(B\) 在图上钦定了一条边 \(i\to p_i\) 必须选,那么在 \(A\) 中选的就应该是:

  • \(p_i\to \cdots\to i\) 的一条链,其中每一条边后者是前者的倍数,且边权为 \(-C\)
  • 发现剩下的点在 \(A\) 中都只能选自环(从大往小推每个点即可),边权为 \(1-C\)

假设 \(f_n\) 表示链尾是 \(n\) 的所有方案的和,可以得到:

\[f_n=C(1-C)^{n-1}+\sum_{d|n\land d\neq n}\frac{C}{1-C}f_d \]

其中前者是链长为 \(1\) 的方案,后者是链长大于 \(1\) 的方案,枚举链尾的前一个数 \(d\),注意环个数的变化所带来的正负号影响。

我们只需要求出 \(\sum_{i=1}^n f_n\) 即可,发现这也是可以用类似杜教筛的方法来做的:

\[\begin{aligned} \sum_{i=1}^nf_i&=\sum_{i=1}^n\left(C(1-C)^{n-1}\sum_{d|i,d\neq i}\frac{C}{1-C}f_d\right)\\ &=nC(1-C)^{n-1}+\sum_{d=1}^nf_d\cdot\frac{C}{1-C}\left(\left\lfloor\frac{n}{d}\right\rfloor-1\right) \end{aligned} \]

有一个问题是如何线性预处理小范围的 \(\sum_{i=1}^n f_i\)。发现 \(f_i\) 只与 \(i\) 的因子的形态有关:具体来说,设 \(i=p_1^{a_1}\cdots p_k^{a_k}\),那么 \(f_i\) 只与可重集 \(\{a_i\}\) 有关。

那么我们为每种 \(\{a_i\}\) 选一个符合它的最小的数作为代表元(具体来说,将 \(a_i\) 从大到小排序后,取代表元为 \(2^{a_1}3^{a_2}\cdots\)),并只对这些代表元暴力求。

代表元数量很小(实测小于等于 \(2\times 10^7\) 的代表元只有不到 \(600\) 个,所以对它们每个都根号暴力求都没问题),这样就能做到线性预处理 \(O(n^{\frac{2}{3}})\) 内的前缀和了。

时间复杂度 \(O(n^{\frac{2}{3}})\)

#include<bits/stdc++.h>

#define ll long long

using namespace std;

namespace modular
{
	const int mod=998244353;
	inline int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
	inline int dec(int x,int y){return x-y<0?x-y+mod:x-y;}
	inline int mul(int x,int y){return 1ll*x*y%mod;}
	inline void Add(int &x,int y){x=x+y>=mod?x+y-mod:x+y;}
	inline void Dec(int &x,int y){x=x-y<0?x-y+mod:x-y;}
	inline void Mul(int &x,int y){x=1ll*x*y%mod;}
}using namespace modular;

inline int poww(int a,ll b)
{
	int ans=1;
	while(b)
	{
		if(b&1) ans=mul(ans,a);
		a=mul(a,a);
		b>>=1;
	}
	return ans;
}

inline int read()
{
	int x=0,f=1;
	char ch=getchar();
	while(ch<'0'||ch>'9')
	{
		if(ch=='-') f=-1;
		ch=getchar();
	}
	while(ch>='0'&&ch<='9')
	{
		x=(x<<1)+(x<<3)+(ch^'0');
		ch=getchar();
	}
	return x*f;
}

const int n2=300000,n3=20000000;
const int N2=n2+10,N3=n3+10;

ll n;
int C,K,coef;
int idx,id1[N2],id2[N2];
int cnt,prime[N3],mup[N3],omega[N3],sign[N3];
bool notprime[N3];
int sumf[N2],pref[N3],ff[N3];

int &id(ll x){return x<=n2?id1[x]:id2[n/x];}

void init()
{
	mup[1]=1,omega[1]=0;
	for(int i=2;i<=n3;i++)
	{
		if(!notprime[i])
		{
			prime[++cnt]=i;
			mup[i]=i,omega[i]=1;
		}
		for(int j=1,v;j<=cnt&&(v=i*prime[j])<=n3;j++)
		{
			notprime[v]=1;
			if(!(i%prime[j]))
			{
				mup[v]=mup[i],omega[v]=omega[i];
				break;
			}
			mup[v]=mup[i]*prime[j],omega[v]=omega[i]+1;
		}
	}
	vector<int> sp;
	sp.push_back(1);
	for(int i=1;i<=cnt;i++)
	{
		sp.push_back(sp.back()*prime[i]);
		if(sp.back()>n3) break;
	}
	sign[1]=1,ff[1]=K;
	for(int i=2;i<=n3;i++)
	{
		sign[i]=sign[i/mup[i]]*sp[omega[i]];
		if(i!=sign[i])
		{
			ff[i]=ff[sign[i]];
			continue;
		}
		ff[i]=K;
		for(int j=1;j*j<=i;j++)
		{
			if(!(i%j))
			{
				Add(ff[i],mul(coef,ff[j]));
				if(j>1&&j*j!=i) Add(ff[i],mul(coef,ff[i/j]));
			}
		}
	} 
	for(int i=1;i<=n3;i++) pref[i]=add(pref[i-1],ff[i]);
}

int solve(ll n)
{
	if(n<=n3) return pref[n];
	if(sumf[id(n)]!=-1) return sumf[id(n)];
	int ans=mul(n%mod,K);
	for(ll l=1,r;;l=r+1)
	{
		r=n/(n/l);
		if(r==n) break;
		Add(ans,mul(mul(coef,(n/l-1)%mod),dec(solve(r),solve(l-1))));
	}
	return sumf[id(n)]=ans;
}

int main()
{
	cin>>n>>C;
	if(C==1)
	{
		puts(n<=2?"1":"0");
		return 0;
	}
	K=mul(C,poww(dec(1,C),n-1));
	coef=mul(C,poww(dec(1,C),mod-2));
	for(ll l=1,r;l<=n;l=r+1)
		r=n/(n/l),id(n/l)=++idx;
	init();
	memset(sumf,-1,sizeof(sumf));
	int res=solve(n);
	printf("%d\n",add(poww(dec(1,C),n),res));
	return 0;
}
/*
2 2
*/
posted @ 2022-10-31 07:47  ez_lcw  阅读(43)  评论(0编辑  收藏  举报