1. 程式人生 > 實用技巧 >盧卡斯Lucas & 擴充套件盧卡斯exLucas

盧卡斯Lucas & 擴充套件盧卡斯exLucas

盧卡斯定理

適用條件:

只算數量較少的C(n, m), 模數P較小(約1e5),且P為質數。

遞推式:C(n, m) = C(n/P, m/P) * C(n%P, m%P) % P;

當 n < m 的時候, \(C_n^m = 0\)

當 n == 0 或者 m == 0 的時候就可以停止了。

Code:

//階乘已經預處理好,逆元用快速冪
long long Lucas(long long n, long long m) {//n > m 
	if (m == 0 || n == 0)	return 1;
	long long res = Lucas(n / p, m / p) * C(n % p, m % p) % p;
	return res;
}

擴充套件盧卡斯

適用範圍:

只求少量組合數,m, n巨大,但模數 \(M = \prod p^q\) ,且 \(p^q\) 不大。

複雜度:\(O(\log_Pn * (p^q))\)

擴充套件盧卡斯實際上還可以搞任何階乘除階乘的類似問題。比如:P2183 [國家集訓隊]禮物。主要思路是對每個 \(p_q\) 單獨考慮,對因子 \(p\) 和其它因子單獨考慮。提出來 \(p\) 的倍數,剩下的利用暴力字首積搞,\(p\) 的倍數部分遞迴搞。

Code:

//P2183禮物
int realP;
int n, m, w[10];
int p[44], c[44], ptot, rp[44];
inline void Div(int x) {
	for (int i = 2; i * i <= x; ++i) {
		if (x % i == 0) {
			p[++ptot] = i; rp[ptot] = 1;
			while (x % i == 0)	x /= i, ++c[ptot], rp[ptot] *= i;
		}
	}
	if (x > 1)	p[++ptot] = x, c[ptot] = 1, rp[ptot] = x;
}

inline int quickpow(int x, int k, int P)

int exgcd(int a, int b, int &x, int &y)

int calc(int a, int b) {//ax + by = 1 (gcd(a,b)=1), return x
	int x, y; exgcd(a, b, x, y);
	return (x % b + b) % b;
}

int nw;
int yu, ct;
int Prod[101000];
int fakeans[44];
void sol(int n) {
	if (!n)	return ;
	yu = yu * quickpow(Prod[rp[nw]], n / rp[nw], rp[nw]) % rp[nw] * Prod[n % rp[nw]] % rp[nw];
	ct += n / p[nw];
	sol(n / p[nw]);
}

inline void merge() {
	int ans = 0;
	for (int i = 1; i <= ptot; ++i) {
		int tmp = calc(realP / rp[i], rp[i]);
		ans = (ans + fakeans[i] * (realP / rp[i]) % realP * tmp) % realP;
	}
	printf("%lld\n", (ans % realP + realP) % realP);
}

signed main() {//calculate (n!)/(w_1! w_2!...)
	Div(realP);
	for (int i = 1; i <= ptot; ++i) {
		nw = i;
		yu = 1, ct = 0;
		Prod[0] = 1;
		for (int j = 1; j <= rp[i]; ++j)
			if (j % p[i])	Prod[j] = Prod[j - 1] * j % rp[i];
			else	Prod[j] = Prod[j - 1];
		sol(n);
		int memo = yu, memoct = ct; yu = 1; ct = 0;
		for (int j = 1; j <= m; ++j)	sol(w[j]);
		yu = calc(yu, rp[nw]) * memo % rp[i];
		fakeans[i] = yu * quickpow(p[i], memoct - ct, rp[i]) % rp[i];
		yu = 1; ct = 0;
	}
	merge();
}