题解 P2480 [SDOI2010]古代猪文

文章介绍了一道涉及数论知识的编程题,要求计算特定指数的模运算结果。解题策略包括应用欧拉定理和其推论简化表达式,利用Lucas定理处理组合数取模,以及通过中国剩余定理解决大数模运算。对于大质数模,首先进行质因数分解,然后分别求解模每个质因子的余数,最后使用中国剩余定理合并结果。代码实现中包含了关键函数如快速幂和Lucas定理计算。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

题意

g∑d∣nCnd mod 999911659g^{\sum\limits_{d | n} C_n^d} \bmod 999911659gdnCndmod999911659

n,g≤109n, g \le 10^9n,g109

一道非常好的数论题,用到了基本所有的基础数论知识。

需要使用到的数论知识

  • 欧拉定理
  • 逆元
  • lucas 定理
  • 中国剩余定理

对式子变形

先令 ∑d∣nCnd=x\sum\limits_{d | n} C_n^d = xdnCnd=x

首先整个式子是可以用快速幂搞出来的,所以我们考虑怎么求出那个巨大的系数。

由于 999911659999911659999911659 是质数,那么考虑使用欧拉定理和它的推论来化简。

欧拉定理:若 gcd⁡(a,m)=1\gcd(a, m) = 1gcd(a,m)=1,则 aφ(m)≡1(modm)a^{\varphi(m)} \equiv 1 \pmod maφ(m)1(modm)

推论:

ab≡{ab mod φ(m),gcd⁡(a,m)=1ab,b<φ(m)ab mod φ(m)+φ(m),b≥φ(m)a^b \equiv \begin{cases} a^{b \bmod \varphi(m)} & , \gcd(a, m) = 1 \\ a^b &, b < \varphi(m) \\ a^{b \bmod \varphi(m) + \varphi(m)} &, b \ge \varphi(m) \end{cases}ababmodφ(m)ababmodφ(m)+φ(m),gcd(a,m)=1,b<φ(m),bφ(m)

那么原式就变成了

gx mod 999911658 mod 999911659g^{x \bmod 999911658} \bmod 999911659gxmod999911658mod999911659

分解模数

求组合数取模考虑使用 lucas 定理,但这个模数太大了,lucas 定理的时间复杂度是 O(plog⁡np)\mathcal{O}(p \log_np)O(plognp) 的,会爆掉。

考虑对其分解质因数。

999911658=2×3×4679×35617999911658 = 2 \times 3 \times 4679 \times 35617999911658=2×3×4679×35617

这个数非常非常的好,因为每个质因数(也就是后面的模数)的指数都是 111,可以直接上中国剩余定理。

如果不是上面这样,还得用扩展中国剩余定理。

我们能用 lucas 分别求出那个系数模上面四个质数的余数 a1∼a4a_1 \sim a_4a1a4

{x≡a1(mod2)x≡a2(mod3)x≡a3(mod4679)x≡a4(mod35617)\begin{cases} x \equiv a_1 \pmod 2 \\ x \equiv a_2 \pmod 3 \\ x \equiv a_3 \pmod {4679} \\ x \equiv a_4 \pmod {35617} \end{cases}xa1(mod2)xa2(mod3)xa3(mod4679)xa4(mod35617)

用 CRT 求出结果

套模板即可。

一个细节

上述算法最开始的化简的要求是 gcd⁡(a,m)=1\gcd(a, m) = 1gcd(a,m)=1

数据范围内有可能出现 g=999911659g = 999911659g=999911659 的情况,需要特判。

代码实现

#include <bits/stdc++.h>
using namespace std;

typedef long long LL;
typedef pair<int, int> PII;

const int inf = 0x3f3f3f3f;
const LL infLL = 0x3f3f3f3f3f3f3f3fLL;

const int P = 999911658;
int mi[] = {2, 3, 4679, 35617};
int a[4];

int n, g;

int qpow(int a, int b, int P)
{
	int res = 1;
	while (b)
	{
		if (b & 1) res = (LL)res * a % P;
		a = (LL)a * a % P;
		b >>= 1;
	}
	return res;
}

int C(int n, int m, int P)
{
	if (n < m) return 0;
	if (m > n - m) m = n - m;
	int c1 = 1, c2 = 1;
	for (int i = 1; i <= m; i ++ )
	{
		c1 = (LL)c1 * (n - i + 1) % P;
		c2 = c2 * i % P;
	}
	return (LL)c1 * qpow(c2, P - 2, P) % P;
}
int lucas(int n, int m, int P)
{
	return !m ? 1 : (LL)C(n % P, m % P, P) * lucas(n / P, m / P, P) % P;
}

int crt()
{
	int m = P;
	int res = 0;
	for (int i = 0; i < 4; i ++ )
	{
		int M = m / mi[i];
		res = (res + (LL)M * qpow(M, mi[i] - 2, mi[i]) * a[i]) % P;
	}
	return res;
}

int main()
{
	cin >> n >> g;
	if (g % (P + 1) == 0) return puts("0") & 1;
	for (int i = 0; i < 4; i ++ )
	{
		for (int j = 1; j * j <= n; j ++ )
		{
			if (n % j) continue;
			a[i] = ((LL)a[i] + lucas(n, j, mi[i]));
			if (j == n / j) continue;
			a[i] = ((LL)a[i] + lucas(n, n / j, mi[i]));
		}
	}
	cout << qpow(g, crt(), P + 1) << endl;
	return 0;
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值