数论,主要用于解决和质数,余数,整数相关的问题。

上面那句话似乎是病句,但不用在意。

本博客不会讲解最基础的逆元求法与裴蜀定理,请读者自学后再来观看。


扩展欧几里得

扩展欧几里得主要用于寻找二元一次不定方程的解。例如以下的式子:

\(ax + by = c\)

其中a, b为常数

首先我们解决这个问题的简化版本。根据裴蜀定理,如果 \(gcd(a,b)|c\) 则原方程有解:

\(ax+by=gcd(a,b)\)

我们发现如果按照欧几里得的步骤对a, b进行转化,最后会达成这样的情况:

\(a’x’ + 0y’ = gcd(a,b)\)

此时根据欧几里得的步骤,\(a’=gcd(a,b)\) 。所以 \(x’=1\)

这启发我们在返回时更新x, y的值。不妨对欧几里得的一个步骤进行推导:

\(ax+by=gcd(a,b)→bx’+(a\;mod\; b)y’=gcd(a,b)\)

\(a\; mod\;b\) 变为 \(a – \lfloor \frac{a}{b} \rfloor * b\),然后可得:

\(bx’ + (a – \lfloor \frac{a}{b} \rfloor * b)y’\)

换一下顺序:

\(ay’ + b(x’ – \lfloor \frac{a}{b} \rfloor * y’)\)

然后可以得出对应关系。之后便可以递归求解。

__int128 exgcd(__int128 &x, __int128 &y, __int128 a, __int128 b)
{
	if(b == 0) 
	{
		x = 1;y = 0;
		return a;	
	}
	__int128 p = exgcd(x, y, b, a % b);
	__int128 x0 = x, y0 = y;
	x = y0, y = x0 - (a / b) * y0;
	return p;
}

注意:尽管在数学上当b=0时y可以取任意数,但为了之后的回溯求解,y只能取0

时间复杂度与欧几里得算法一样,为log级别


欧拉函数

主要用于优化连快速幂也解决不了的幂运算

介绍

欧拉函数用 \(φ\) 表示。\(φ(n)\) 表示不大于 \(n\) 的数中与 \(n\) 互质的正整数个数。

欧拉函数是一个积性函数,所以她有性质:当 \(gcd(a,b)=1\) 时, \(φ(ab)=φ(a)φ(b)\)

性质

性质一:\(φ(1) = 1\)

性质二:当 \(n\) 为质数时,\(φ(n)=n-1\)

以上两条可以通过性质得出

性质三:当 \(p\) 为质数时, \(φ(p^k)=p^k-p^{k-1}=p^{k-1}(p-1)\) 。实际上,\(p^k\) 的因数只有 \(p^{k-1}\) 个减去就可以了。

性质四:\(φ(n)= n * \prod_{i=1}^s{\dfrac{p_i-1}{p_i}}\), 其中数组 \(p\)\(n\) 的质因子

证明如下:

根据积性函数性质有:

\(φ(n)= \prod_{i=1}^s{φ(p_i^{k_i})}\),其中 \(n= \prod_{i=1}^s{p_i^{k_i}}\)

根据性质三:

\(φ(n)= \prod_{i=1}^s{p_i^{k_i-1}(p_i-1)}\)

发现 \(\prod_{i=1}^s{p_i^{k_i-1}}\)\(n\) 只有一个 \(\prod_{i=1}^s{p_i}\) 的差距。所以可以化为:

\(φ(n)= (n * \prod_{i=1}^s{\dfrac{1}{{p_i}}}) \prod_{i=1}^s{{p_i-1}}\)

所以有

\(φ(n)= n * \prod_{i=1}^s{\dfrac{p_i-1}{p_i}}\)

很简单吧?

求法

对于单个数的欧拉函数,我们可以直接通过性质四求解。可以用 Pollard Rho 优化。时间复杂度大约为 log 级别。

int ol(int p)
{
	int ans = p;
	for(int i = 2;i * i <= p;i ++)
	{
		if(p % i == 0) ans = ans / i * (i - 1);
		while(p % i == 0) p /= i;
	}
	if(p > 1) ans = ans / p * (p - 1);
	return ans;
}

如果我们想线性求出从 1 开始的自然数的欧拉函数,就要依靠欧拉筛了原汤化原食

在欧拉筛中,每一个合数会在她最小的质因子处被标记。这启发我们将 \(n\) 拆成最小的质因子 \(p\)\(n’=\dfrac{n}{p}\)

因为 \(p\) 是质数,所以只有 \(n’\)\(p\) 的倍数,和 \(n’\)\(p\) 互质两种情况。

先考虑第一种:当 \(n’\)\(p\) 倍数时,就代表 \(n’\)\(n\) 质因子完全一样。于是就有下列转化:

\(φ(n)= n * \prod_{i=1}^s{\dfrac{p_i-1}{p_i}} = p * n’ * \prod_{i=1}^s{\dfrac{p_i-1}{p_i}} = p * φ(n’)\)

第二种情况:当 \(n’\)\(p\) 互质时,运用积性函数的性质可以直接得出:\(φ(n)=φ(p)φ(n’)\)

于是可以得出线性求法。

int ol[100005], zhi[100005], cnt;
void OL(int n)
{
	ol[1] = 1;
	for(int i = 2;i <= n;i ++)
	{
		//当走到 i 时,若 i 为合数,则 ol[i] 已被求出,要做的只是更新。
		if(!ol[i])
		{
			ol[i] = i - 1;
			zhi[++ cnt] = i;
		} 
		for(int _ = 1;_ <= cnt && i * zhi[_] <= 100000;_ ++)
		{
			int p = zhi[_];
			if(i % zhi[_]) ol[i * zhi[_]] = ol[i] * (zhi[_] - 1);
			else
			{
				ol[i * zhi[_]] = ol[i] * zhi[_];
			}
		}
	}
}

欧拉定理

欧拉定理原版:若 \(gcd(a,m) = 1\)\(a^{φ(m)} \equiv 1(mod \;m)\)

欧拉定理扩展:若 \(b > φ(m)\) ,则 $a^b \equiv a^{(b;mod; m) + m} (mod ; m) $

实际上扩展欧拉定理用的更多,因为大多数时候模数并不大。

证明之后再补。

例题

[SDOI2010]古代猪文

题目描述

猪王国的文明源远流长,博大精深。

iPig 在大肥猪学校图书馆中查阅资料,得知远古时期猪文文字总个数为 \(n\)。当然,一种语言如果字数很多,字典也相应会很大。当时的猪王国国王考虑到如果修一本字典,规模有可能远远超过康熙字典,花费的猪力、物力将难以估量。故考虑再三没有进行这一项劳猪伤财之举。当然,猪王国的文字后来随着历史变迁逐渐进行了简化,去掉了一些不常用的字。

iPig 打算研究古时某个朝代的猪文文字。根据相关文献记载,那个朝代流传的猪文文字恰好为远古时期的 \(1/k\),其中 \(k\)\(n\) 的一个正约数(可以是 \(1\)\(n\))。不过具体是哪 \(1/k\),以及 \(k\) 是多少,由于历史过于久远,已经无从考证了。

iPig 觉得只要符合文献,每一种 \(k|n\) 都是有可能的。他打算考虑到所有可能的 \(k\)。显然当 \(k\) 等于某个定值时,该朝的猪文文字个数为 \(n/k\)。然而从 \(n\) 个文字中保留下 \(n/k\) 个的情况也是相当多的。iPig 预计,如果所有可能的 \(k\) 的所有情况数加起来为 \(p\) 的话,那么他研究古代文字的代价将会是 \(g^p\)

现在他想知道猪王国研究古代文字的代价是多少。由于 iPig 觉得这个数字可能是天文数字,所以你只需要告诉他答案除以 \(999911659\) 的余数就可以了。

提示

数据规模与约定

  • 对于 \(10\%\) 的数据,\(1\le n \le 50\)
  • 对于 \(20\%\) 的数据,\(1\le n \le 1000\)
  • 对于 \(40\%\) 的数据,\(1\le n \le 10^5\)
  • 对于 \(100\%\) 的数据,\(1\le n,g \le 10^9\)

要我们求的式子实际上为:\(g^{\sum_{d|n}\binom{n}{d} }\)

这个式子可以用欧拉定理优化为:

\(g^{\sum_{d|n}\binom{n}{d} mod 999911658}\)

实际我们要求的就是

\(\sum_{d|n}\binom{n}{d} \;mod\;999911658\)

之后交给快速幂就好。

由于有 C ,想到 \(lucas\) 定理。但这个东西直接用会挂。于是考虑缩减模数。

发现 999911658 的质因子很少,只有2,3,4679,35617。所以我们可以将模数为这四个数时的结果求出来,再用中国剩余定理合并即可。

#include <bits/stdc++.h>
#define ll long long
using namespace std;
ll n, g;
ll a[4] = {2, 3, 4679, 35617}, b[4], jie[500005];
ll qpow(ll u, ll v, ll p)
{
	ll ans = 1;
	while(v)
	{
		if(v & 1) ans *= u, ans %= p;
		u *= u;u %= p;
		v >>= 1;
	}
	return ans;
}
ll C(ll n, ll m, ll p)
{
	if(n < m) return 0;
	return (((jie[n] * qpow(jie[m], p - 2, p)) % p) * qpow(jie[n - m], p - 2, p)) % p;
}
ll lucas(ll n, ll m, ll p)
{
	if(n < m) return 0;
	if(n == 0) return 1;
	return (C(n % p, m % p, p) * lucas(n / p, m / p, p)) % p;
}
ll exgcd(ll &x, ll &y, ll a, ll b)
{
	if(b == 0)
	{
		x = 1, y = 0;
		return a;
	}
	ll gcd = exgcd(x, y, b, a % b);
	ll x0 = x, y0 = y;
	x = y0; y = x0 - (a / b * y0);
	return gcd;
}
ll excrt(ll a[], ll b[], ll n)
{
	ll ans = b[0], lcm = a[0];
	for(int i = 1;i < n;i ++)
	{
		b[i] = ((b[i] - ans) % a[i] + a[i]) % a[i];
		ll x, y;
		ll d = exgcd(x, y, lcm, a[i]);
		ans = (ans + x * b[i] * (lcm / d));
		lcm = (lcm / d * a[i]);
		ans = ((ans % lcm) + lcm) % lcm;
	}
	return ans;
}
ll query(ll p)
{
	ll ans = 0;
	for(ll i = 1;i * i <= n;i ++)
	{
		if(n % i) continue;
		ans = (ans + lucas(n, i, p)) % p;
		if(i * i < n)  ans = (ans + lucas(n, n / i, p)) % p;
	}
	return ans;
}
void init(ll p)
{
	jie[0] = jie[1] = 1;
	for(int i = 2;i <= p;i ++) jie[i] = (jie[i - 1] * i) % p;
}
int main()
{
	cin >> n >> g;
	if(g % 999911659 == 0)
	{
		cout << 0 << endl;
		return 0;
	}
	ll p = 0;
	for(int i = 0;i < 4;i ++)
	{
		init(a[i]);
		b[i] = query(a[i]);
	}
	ll P = excrt(a, b, 4);
	cout << qpow(g, P, 999911659) << endl; 
	return 0;
}

扩展中国剩余定理

由于扩展中国剩余定理与中国剩余定理之间做法几乎无关联,且两者复杂度一样,但扩展中国剩余定理适用性更强,所以此处只介绍扩展中国剩余定理。

扩展中国剩余定理用于解多个同余方程。例如

\[\begin{cases} x \equiv b_1\ ({\rm mod}\ a_1) \\ x\equiv b_2\ ({\rm mod}\ a_2) \\ … \\ x \equiv b_n\ ({\rm mod}\ a_n)\end{cases} \]

设我们当前求出了满足前 i – 1 个方程的值 x ,要据此求出满足前 i 个方程的值 x。

\(m = gcd(a_1,a_2…,a_{i-1})\),则因为要满足前i-1个方程,所以最终解必须为 \(x+tm\) 的形式(t为常数)

那么实际上方程为:\(x+tm \equiv b_i(mod\;a_i)\)

既然x为定值,那么可以把x移动到右边,使方程变为\(tm \equiv b_i-x(mod\;a_i)\)

于是可以化为:\(tm+pa_i=b_i-x\)。即一个以t,p为未知数的二元一次不定方程。用上文提到的扩展欧几里得求解即可。

当然有可能出现求解过程中方程无解的情况。此时方程组当然也无解。

时间复杂度为O(nlogn)。其中第二个n其实是值域

long long exchina()
{
	__int128 ans = b[1], lcm = a[1];
	for(int i = 2;i <= n;i ++)
	{
		b[i] = ((b[i] - ans) % a[i] + a[i]) % a[i];//防止因b[i]<ans导致的负数情况
		__int128 x, y;
		__int128 d = exgcd(x, y, lcm, a[i]);
		if(b[i] % d) return -1;
		ans = (ans + x * b[i] * lcm / d);
		lcm = lcm * (a[i] / d);
		ans = ((ans % lcm) + lcm) % lcm;
	}
	return ans;
}

原文地址:http://www.cnblogs.com/closureshop/p/16866437.html

1. 本站所有资源来源于用户上传和网络,如有侵权请邮件联系站长! 2. 分享目的仅供大家学习和交流,请务用于商业用途! 3. 如果你也有好源码或者教程,可以到用户中心发布,分享有积分奖励和额外收入! 4. 本站提供的源码、模板、插件等等其他资源,都不包含技术服务请大家谅解! 5. 如有链接无法下载、失效或广告,请联系管理员处理! 6. 本站资源售价只是赞助,收取费用仅维持本站的日常运营所需! 7. 如遇到加密压缩包,默认解压密码为"gltf",如遇到无法解压的请联系管理员! 8. 因为资源和程序源码均为可复制品,所以不支持任何理由的退款兑现,请斟酌后支付下载 声明:如果标题没有注明"已测试"或者"测试可用"等字样的资源源码均未经过站长测试.特别注意没有标注的源码不保证任何可用性