hdu 3439 Sequence 错排公式 + 排列组合+CRT %不是质数

题意:

1~n的排列有n!种，定义D(某个排列)表示为排列中不动点的个数。例如D({1,2,3})=3,D({1,3,2})=1。问1~n的排列中，不动点个数为k的有多少个，对m取余 (T<=500 0<=k<=n<=1e9 1<=m<=1e5, n!=0)

思路:

其实就是一个组合数加错排？ C(n,k)*F(n-k) 选出k个不动,剩下的错排公式.

关键是取%不是个质数了.还是老方法,对m分解质因数分别求然后CRT合成.

观察发现n这么大错排怎么错排呢？发现模的m很小啊,这种情况下一般有个通用的解法就是打表找规律.可以发现2*m为一个循环节.

另外附上错排公式的两种求法:

 F[N] = (N-1)*(F[N-1]+F[N-2])

F[N] = N*F[N-1]+(-1)^N

#include<bits/stdc++.h>using namespace std;typedef long long ll;const int mod=1e9+7;const int maxn = 5e6 + 10;ll p[233],pc[233],M[233],a[233];ll num,fac[maxn];ll n,m,k;ll qmod(ll a,ll b,ll mod){ll res = 1;while(b){if(b&1)res = res * a % mod;a = a * a % mod;b >>= 1;}return res;}void exgcd(ll a,ll b,ll &x,ll &y){    if(!b){x=1;y=0;return;}    else {exgcd(b,a%b,y,x);y-=x*(a/b);}}ll inv(ll x,ll p)//扩展欧几里得求逆元 {    ll xx,yy;    exgcd(x,p,xx,yy);    return ((xx%p)+p)%p;}ll calc(ll x,ll p,ll s){if(x < p) return fac[x];num += x/p;return fac[x % s]*qmod(fac[s-1],x/s,s)%s*calc(x/p,p,s)%s;}ll work(ll p,ll s)//计算 C(n,n-k) {fac[0] = 1;for(ll i = 1;i < s ;++i)fac[i] = fac[i-1]*(i%p==0?1:i)%s;num = 0;ll n1 = calc(n,p,s);// n!ll tmp = num;for(ll i = 1;i < s ;++i)//这里的n可能会大于%的m,所以加不加上界都可以 fac[i] = fac[i-1]*(i%p==0?1:inv(i,s))%s;num = 0;ll n2 = calc(k,p,s)*calc(n-k,p,s)%s; n!*(n-k)!return n1 * n2 % s * qmod(p,tmp - num,s) % s;}//错位排列的周期性//打表可得,错排公式,在%mod情况下,每2*mod 为一个周期 ll F(ll x) {ll res = 0;if(x == 0) return 1;x = x % (2*m);if(x == 0) x = 2*m;for(int i = 2;i <= x;++i)res = (res*i + (i % 2 == 0 ?1 : -1))%m;return (res + m ) % m;}int main(){int _;cin>>_;int ca = 1;while(_--){scanf("%lld %lld %lld",&n,&k,&m);ll tmp = m;int len = 0;//分解 m 为 p1^c1*p2^c2..... for(ll i = 2;i*i<=tmp&&tmp>1;++i){if(tmp % i ==0){p[++len] = i;pc[len] = 1;}while(tmp % i == 0){pc[len] *= i;tmp /= i;}}if(tmp > 1){p[++len] = tmp;pc[len] = tmp;}ll ans = 0;for(ll i=1;i<=len;++i)//中国剩余定理合成     {        M[i]=m/pc[i];        a[i]=work(p[i],pc[i]);        ans=(ans+a[i]*M[i]%m*inv(M[i],pc[i])%m)%m;    }    ans = ans * F(n-k) %m;    printf("Case %d: %lld\n",ca++,ans);}return 0;}