自动搬运

查看原文

来自洛谷，原作者为

	orangejuice STATE=HLT

搬运于2025-08-24 22:19:13，当前版本为作者最后更新于2020-09-27 20:45:41，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

---

「USACO 2020 US Open Platinum」Exercise

做法与模数是否是质数无关

问题可能比较复杂，需要多步分析

1.对于一个已知的排列

显然这样的置换会构成若干个环，设每个环长度为$a_i,i\in [1,m]$，显然答案就是$lcm(a_i)$

2.对于已知的$a_i$序列(注意这里说的是有序的)，计算其方案数

考虑已经排列的个数为$i$，加入一个环大小为$j$

为了避免重复，应当固定这个环的初始位置为1号点，其余位置按照原先顺序插入

则方案数可以分为两部分考虑:

2-1.环内排列，固定的环首不可排列，即$(j-1)!$

2-2.剩下的$j-1$一个点位置未知，从未固定的$i+j-1$个点中选择

即$C(i+j-1,j-1)$

所以就是$C(i+j-1,j-1)(j-1)!=\frac{(i+j-1)!}{i!}$

归纳一下，发现更形象的描述就是$\begin{aligned}\frac{n!}{\prod_{i=1}^{m-1} (\sum_{j=1}^i a_j)}\end{aligned}$

也就是每次除掉转移时的大小，将$n!$分成若干段，这似乎有利于理解下面的dp优化

3.计算$lcm$之积

考虑对于每个质因数计算其出现的幂次，注意这个幂次是对于$\varphi$取模的

原先是求恰好包含$x^k$的方案数，得到的$dp$不好优化，考虑转换为:

求质因数$x$出现在答案里的幂次$\ge k$的方案数$F_k$，答案就是$x^{\sum F_k}$

Solution 1

那么反向求解，令$dp_i$表示当前已经确定了$i$个点，没有出现$x^k$倍数大小的联通块

暴力转移，枚举$i$从所有$x^k\not|j$转移过来即可，单次求解复杂度为$O(n^2)$，不可行

优化1:

考虑分解系数$\frac{(i+j-1)!}{i!}$，累前缀和，对于$j$为$x^k$倍数的情况枚举减掉

这样单次求解复杂度为$O(\frac{n^2}{x^k})$，总复杂度为$O(n^2\ln n)$，且不好处理阶乘逆元

优化2:

不枚举$x^i$的倍数，直接再用一个前缀和数组$s_i$记录下来，让$s_i$从$s_{i-x^k}$转移过来即可

如何将系数$\frac{(i+j-1)!}{i!}$分解？

每次$i+j$增大1，就多乘上一个$i+j$即可

当$s_{i}$从$s_{i-x^k}$转移过来时，需要补上$\begin{aligned}\prod_{j=i-x^k+1}^{i-1}j\end{aligned}$

也就是模拟了上面提到的把$n!$分段的过程

这样就去掉了阶乘逆元的求解

Tips:发现需要预处理$T_{i,j}=\prod_{k=i}^j k$，可以滚动一下会快一点，内存为$O(n)$

$$\$$

不同的$x^i$上限为$O(n)$种，实际大概可能是$O(\pi(n)\log \log n=\frac{n\log \log n}{\log n})$?

因此复杂度为$O(n^2)$

可以看到代码还是很简单的

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define Mod2(x) ((x<0)&&(x+=P2))
#define rep(i,a,b) for(int i=a,i##end=b;i<=i##end;++i)
#define drep(i,a,b) for(int i=a,i##end=b;i>=i##end;--i)
 
const int N=7510;
 
int n,P,P2;
int mk[N];
int s[N],T[N],dp[N];
 
ll qpow(ll x,ll k=P-2) {
    ll res=1;
    for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
    return res;
}
 
 
int main(){
    scanf("%d%d",&n,&P),P2=P-1;
    int ans=1;
    rep(i,2,n) if(!mk[i]) {
        for(int j=i+i;j<=n;j+=i) mk[j]=1;
        for(int x=i;x<=n;x*=i) mk[x]=i;
    }
    rep(i,1,n) T[i]=1;
    int r=1;
    rep(i,1,n) r=1ll*r*i%P2;
    rep(x,2,n) {
        rep(j,1,n-x+1) T[j]=1ll*T[j]*(j+x-2)%P2;
        // 滚动求解区间乘积
        if(mk[x]<=1) continue;
        dp[0]=1,s[0]=1;
        int sum=1;
 
        rep(i,1,n) {
            s[i]=0;
            if(i>=x) s[i]=1ll*s[i-x]*T[i-x+1]%P2;
            dp[i]=sum-s[i]; Mod2(dp[i]);
            s[i]=(1ll*s[i]*i+dp[i])%P2;
            sum=(1ll*sum*i+dp[i])%P2;
        }
        ans=1ll*ans*qpow(mk[x],P2+r-dp[n])%P;
    }
    printf("%d\n",ans);
}

Solution 2

为了便于表达，设满足条件为至少出现一个$x$的倍数

实际用min-max容斥确实比较好理解，设对于集合$S$，求其最大值

$\begin{aligned} \max \lbrace S\rbrace =\sum_{T\sube S} (-1)^{|T|+1}\min\lbrace T \rbrace\end{aligned}$

简要证明的话:

把$S$中的元素倒序排成一排分别为$S_i$

对于$S_1$即最大值，显然被计算一次

对于剩下的值$S_i(i>1)$，则它作为最小值产生贡献意味着选的数都在$1-i$内，显然有$2^{i-1}$次为奇数集合大小，$2^{i-1}$为偶数集合大小，两部分抵消

$$\$$

要计算最大值为1的方案数，那么就要计算最小值为1的子集方案数

考虑强制一个子集中每一个环大小均为$x$的倍数，设选出了$i$个这样的环，总大小为$j$的方案数为$dp_{i,j}$

则实际对答案的贡献还要考虑这样的子集出现的次数

考虑选择子集的位置，以及剩下的$n-j$个点任意排布，方案数应该为$C(n,j)\cdot (n-j)!=\frac{n!}{j!}$

如果真的用$dp_{i,j}$，复杂度显然太高，考虑$i$这一维的影响只在于系数$(-1)^{|T|+1}$，可以直接在转移过程中解决

因此可以直接记录大小$i$，从前面转移过来

(可以看到依然需要访问上面提到的$T_{i,j}$，要滚动的话还会更难处理)

这样的$dp$状态有$\frac{n}{x}$种，转移为$(\frac{n}{x})^2$，最后统计复杂度为$O(\frac{n}{x})$

实际上这样的复杂度已经足够了，是$O(\sum \frac{n^2}{i^2}\approx n^2)$

以下是$O(n^2)$内存的代码

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for(int i=a;i<=b;++i)
#define drep(i,a,b) for(int i=a;i>=b;--i)
typedef long long ll;

const int N=7510;

int n,P,P2,T[N][N],mk[N],dp[N];
ll qpow(ll x,ll k) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}
int Calc(int x){
	int m=n/x,s=P2-1;
	rep(i,1,m) {
		s=1ll*s*T[(i-1)*x+1][i*x-1]%P2;
		dp[i]=P2-s;
		s=(1ll*s*i*x+dp[i])%P2;
	}
	int ans=0;
	rep(i,1,m) ans=(ans+1ll*dp[i]*T[i*x+1][n])%P2;
	return ans;
}

int main(){
	scanf("%d%d",&n,&P),P2=P-1;
	rep(i,2,n) if(!mk[i]) {
		for(int j=i;j<=n;j+=i) mk[j]=1;
		for(int j=i;j<=n;j*=i) mk[j]=i;
	}
	rep(i,1,n+1){
		T[i][i-1]=1;
		rep(j,i,n) T[i][j]=1ll*T[i][j-1]*j%P2;
	}
	int ans=1;
	rep(i,2,n) if(mk[i]>1) ans=ans*qpow(mk[i],Calc(i))%P;
	printf("%d\n",ans);
}

$$\$$

最终优化

可以看到，这个做法和Sol1的转移有十分的相同之处，因此考虑用同样的方法优化掉

但是预处理系数的部分依然是$O(n^2)$的，如何解决呢?

1.线段树大法

预处理复杂度为$O(n)$，查询复杂度为$O(\log n)$，总复杂度$O(n\log^2 n)$

空间复杂度为$O(n)$

Oh这个代码是ZKW线段树

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for(int i=a;i<=b;++i)
typedef long long ll;

const int N=7510;

int n,P,P2;
struct Tree {
	int s[N<<2],bit;
	void Build() {
		for(bit=1;bit<=n;bit<<=1);
		rep(i,1,n) s[i+bit]=i;
		for(int i=bit;i>=1;--i) s[i]=1ll*s[i<<1]*s[i<<1|1]%P2;
	}
	int Que(int l,int r){
		if(l>r) return 1;
		if(l==r) return l;
		int res=1;
		for(l+=bit-1,r+=bit+1;l^r^1;l>>=1,r>>=1){
			if(~l&1) res=1ll*res*s[l^1]%P2;
			if(r&1) res=1ll*res*s[r^1]%P2;
		}
		return res;
	}
} T;

int mk[N];
ll qpow(ll x,ll k) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}

int dp[N];
int Calc(int x){
	int m=n/x,s=P2-1;
	rep(i,1,m) {
		s=1ll*s*T.Que((i-1)*x+1,i*x-1)%P2;
		dp[i]=P2-s;
		s=(1ll*s*i*x+dp[i])%P2;
	}
	int ans=0;
	rep(i,1,m) ans=(ans+1ll*dp[i]*T.Que(i*x+1,n))%P2;
	return ans;
}

int main(){
	scanf("%d%d",&n,&P),P2=P-1;
	T.Build();
	rep(i,2,n) if(!mk[i]) {
		for(int j=i;j<=n;j+=i) mk[j]=1;
		for(int j=i;j<=n;j*=i) mk[j]=i;
	}
	int ans=1;
	rep(i,2,n) if(mk[i]>1) ans=ans*qpow(mk[i],Calc(i))%P;
	printf("%d\n",ans);
}

2.模数分解+CRT(Chinese Remainder Theory)

分解后，可以用模逆元处理，然后就直接做，最后CRT合并一下，其实我并不会实现。。。

3.猫树~~(嘿嘿)~~

这是一个$O(n\log n)$预处理，$O(1)$查询的数据结构，空间复杂度为$O(n\log n)$

因此复杂度为$O(n\log n)$

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for(int i=a;i<=b;++i)
#define drep(i,a,b) for(int i=a;i>=b;--i)
typedef long long ll;

const int N=7510;

int n,P,P2;
struct SuckCat{
	int s[14][8200],Log[8200],bit;
	void Build() {
		for(bit=1;bit<=n;bit<<=1);
		rep(i,2,bit) Log[i]=Log[i>>1]+1;
		for(int l=1,d=0;l*2<=bit;l<<=1,d++){
			for(int i=l;i<=bit;i+=l*2){
				s[d][i-1]=i-1,s[d][i]=i;
				drep(j,i-2,i-l) s[d][j]=1ll*s[d][j+1]*j%P2;
				rep(j,i+1,i+l-1) s[d][j]=1ll*s[d][j-1]*j%P2;
			}
		}
	}
	int Que(int l,int r){
		if(l>r) return 1;
		if(l==r) return l;
		int d=Log[l^r];
		return 1ll*s[d][l]*s[d][r]%P2;
	}
} T;

int mk[N];
ll qpow(ll x,ll k) {
	ll res=1;
	for(;k;k>>=1,x=x*x%P) if(k&1) res=res*x%P;
	return res;
}

int dp[N];
int Calc(int x){
	int m=n/x,s=P2-1;
	rep(i,1,m) {
		s=1ll*s*T.Que((i-1)*x+1,i*x-1)%P2;
		dp[i]=P2-s;
		s=(1ll*s*i*x+dp[i])%P2;
	}
	int ans=0;
	rep(i,1,m) ans=(ans+1ll*dp[i]*T.Que(i*x+1,n))%P2;
	return ans;
}

int main(){
	scanf("%d%d",&n,&P),P2=P-1;
	T.Build();
	rep(i,2,n) if(!mk[i]) {
		for(int j=i;j<=n;j+=i) mk[j]=1;
		for(int j=i;j<=n;j*=i) mk[j]=i;
	}
	int ans=1;
	rep(i,2,n) if(mk[i]>1) ans=ans*qpow(mk[i],Calc(i))%P;
	printf("%d\n",ans);
}