【LOJ #3284】「USACO 2020 US Open Platinum」Exercise(容斥 / 组合数学 / DP)

传送门

首先 n n n个数的置换个数是 n ! n! n!
证明可以用生成函数
循环的 e g f   f = ∑ i ( i − 1 ) ! i ! x i = ∑ i x i i egf\ f=\sum_i\frac{(i-1)!}{i!}x^i=\sum_i\frac{x^i}{i} egf f=ii!(i1)!xi=iixi
ln ⁡ ( 1 1 − x ) = ∑ i x i i = f \ln(\frac{1}{1-x})=\sum_i\frac{x^i}{i}=f ln(1x1)=iixi=f
exp ⁡ f = 1 1 − x = ∑ i x i \exp f=\frac{1}{1-x}=\sum_{i}x^i expf=1x1=ixi
因为是 e g f egf egf乘上 i ! i! i!即可


首先可以枚举每个循环的大小 a a a
那么就是这样一个东西
∏ ∑ a = n ( l c m ( a ) ) n ! ∏ a ! \prod_{\sum_a=n}(\mathrm{lcm}(a))^{\frac{n!}{\prod a!}} a=n(lcm(a))a!n!

然后考虑把 l c m lcm lcm拆成 ∏ i p i k i \prod_{i}p_i^{k_i} ipiki
然后考虑每个质数次幂的贡献
∏ j p j ∑ k f ( p j , k ) \prod_{j}p_j^{\sum_kf(p_j,k)} jpjkf(pj,k)
= ∏ j p j ∑ k g ( p j k ) =\prod_{j}p_j^{\sum_kg(p_j^k)} =jpjkg(pjk)
f ( p j , k ) f(p_j,k) f(pj,k)表示有多少种 a a a序列满足其中 p j p_j pj最高次幂为 k k k的方案数
g ( p j k ) g(p_j^k) g(pjk)表示有多少种满足存在 p j k p_j^k pjk的倍数
v = p j k v=p_j^k v=pjk
再考虑容斥
f i f_i fi表示所有循环为 v v v的倍数,大小总和为 v i vi vi的方案数
可以枚举当前的循环大小,钦定这个循环第一个点
就是 f i = ∑ j = 1 i ( i v − 1 j v − 1 ) ( j v − 1 ) ! f i − j , f 0 = 1 f_i=\sum_{j=1}^i{iv-1\choose jv-1}(jv-1)!f_{i-j},f_0=1 fi=j=1i(jv1iv1)(jv1)!fij,f0=1
那么枚举 v v v的倍数的大小,剩下随便选
就是 ∑ i ( n i v ) f i ( n − i v ) ! \sum_i{n\choose iv}f_i(n-iv)! i(ivn)fi(niv)!
但是会算重,考虑容斥系数是循环个数
于是把容斥系数带入 f f f
f i = − ∑ j = 1 i ( i v − 1 j v − 1 ) ( j v − 1 ) ! f i − j , f 0 = − 1 f_i=-\sum_{j=1}^i{iv-1\choose jv-1}(jv-1)!f_{i-j},f_0=-1 fi=j=1i(jv1iv1)(jv1)!fij,f0=1
复杂度为 ∑ p i ∑ k ⌊ n p i k ⌋ 2 = n 2 ∑ p i ∑ k ⌊ 1 p i k ⌋ = n 2 π 2 6 \sum_{p_i}\sum_k\lfloor\frac{n}{p_i^k}\rfloor^2=n^2\sum_{p_i}\sum_k\lfloor\frac{1}{p_i^k}\rfloor=n^2\frac{\pi^2}{6} pikpikn2=n2pikpik1=n26π2

#include<bits/stdc++.h>
using namespace std;
#define cs const
#define re register
#define pb push_back
#define pii pair<int,int>
#define ll long long
#define y1 shinkle
#define fi first
#define se second
#define bg begin
cs int RLEN=1<<20|1;
inline char gc(){
    static char ibuf[RLEN],*ib,*ob;
    (ib==ob)&&(ob=(ib=ibuf)+fread(ibuf,1,RLEN,stdin));
    return (ib==ob)?EOF:*ib++;
}
inline int read(){
    char ch=gc();
    int res=0;bool f=1;
    while(!isdigit(ch))f^=ch=='-',ch=gc();
    while(isdigit(ch))res=(res+(res<<2)<<1)+(ch^48),ch=gc();
    return f?res:-res;
}
inline ll readll(){
    char ch=gc();
    ll res=0;bool f=1;
    while(!isdigit(ch))f^=ch=='-',ch=gc();
    while(isdigit(ch))res=(res+(res<<2)<<1)+(ch^48),ch=gc();
    return f?res:-res;
}
inline int readstring(char *s){
	int top=0;char ch=gc();
	while(isspace(ch))ch=gc();
	while(!isspace(ch)&&ch!=EOF)s[++top]=ch,ch=gc();
	s[top+1]='\0';return top;
}
template<typename tp>inline void chemx(tp &a,tp b){a=max(a,b);}
template<typename tp>inline void chemn(tp &a,tp b){a=min(a,b);}
int mod;
inline int add(int a,int b){return (a+b)>=mod?(a+b-mod):(a+b);}
inline int dec(int a,int b){return (a<b)?(a-b+mod):(a-b);}
inline int mul(int a,int b){static ll r;r=(ll)a*b;return (r>=mod)?(r%mod):r;}
inline void Add(int &a,int b){a=(a+b)>=mod?(a+b-mod):(a+b);}
inline void Dec(int &a,int b){a=(a<b)?(a-b+mod):(a-b);}
inline void Mul(int &a,int b){static ll r;r=(ll)a*b;a=(r>=mod)?(r%mod):r;}
inline int ksm(int a,int b,int res=1){for(;b;b>>=1,Mul(a,a))(b&1)&&(Mul(res,a),1);return res;}
inline int Inv(int x){return ksm(x,mod-2);}
inline int fix(ll x){x%=mod;return (x<0)?x+mod:x;}
cs int N=7505;
char xxx;
int c[N][N],n;
int fac[N],f[N],g[N];
int pr[N],tot;
bool vis[N];
char yyy;
inline void init_sieve(int n){
	for(int i=2;i<=n;i++){
		if(!vis[i])pr[++tot]=i;
		for(int j=1;j<=tot&&i*pr[j]<=n;j++){
			vis[i*pr[j]]=1;
			if(i%pr[j]==0)break;
		}
	}
}
inline void init_c(){
	for(int i=0;i<=n;i++){
		c[i][i]=c[i][0]=1;
		for(int j=1;j<i;j++)
		c[i][j]=add(c[i-1][j],c[i-1][j-1]);
	}
	fac[0]=1;
	for(int i=1;i<=n;i++)fac[i]=mul(fac[i-1],i);
}
inline int calc(int v){
	int ln=n/v;
	f[0]=mod-1;
	for(int i=1;i<=ln;i++){
		int res=0;
		for(int j=1;j<=i;j++)
		Dec(res,mul(c[i*v-1][j*v-1],mul(fac[j*v-1],f[i-j])));
		f[i]=res;
	}int res=0;
	for(int i=1;i<=ln;i++){
		Add(res,mul(f[i],mul(fac[n-i*v],c[n][i*v])));
	}return res;
}
int main(){
	#ifdef Stargazer
	freopen("lx.in","r",stdin);
	#endif
	n=read(),mod=read();
	mod=mod-1;
	init_c(),init_sieve(n);
	mod=mod+1;
	int res=1;
	for(int i=1;i<=tot;i++){
		int pw=0;mod=mod-1;
		for(int mt=pr[i];mt<=n;mt*=pr[i])
			Add(pw,calc(mt));
		//,cout<<pw<<'\n';
		mod=mod+1,Mul(res,ksm(pr[i],pw));
	}
	cout<<res<<'\n';
	return 0;
}
  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值