信息学中的数论(一)
做oi题目的时候,遇到数论题会令我兴奋不已。
这一篇让我来聊一聊我学过的gcd,lcm,扩展欧几里得算法,逆元,组合数等。
这篇贴的代码都是未经过编译运行的,所以如果有错或有疑问请评论。
恩
那么什么是数论
和数学有关的非几何都是数论?
嘛,我也不知道定义,那么就草率地认为所有和数学有关的非计算几何知识都是数论吧。
我们先来聊一聊gcd。
这个东西,非常的有用。
它的名字叫最大公约数。
正常人都知道,有一个方法叫辗转相除法(证明略):
int gcd(int a,int b)
{
if(!b)return a;
return gcd(b,a%b);
}
我们默认gcd的时间复杂度是log的,你也可以认为是O(1)的,只不过它的最坏时间复杂度是log的。
不过像我这种渣渣貌似不是很会证明,感受一下是这样的。。
那么gcd讲完了。
我们来看lcm。
lcm的名字叫最小公倍数.
正常人都知道,gcd(a,b)*lcm(a,b)=a*b
那么只要求出gcd,就可以求出lcm了。
我们来证明一下,
假如a和b的gcd是g,那么设a=A*g,b=B*g,易证gcd(A,B)=1,
所以lcm(a,b)=lcm(A*g,B*g)=A*B*g=a*b/g
所以gcd(a,b)*lcm(a,b)=a*b。
这个东西,或许有用。
我们来讲一下扩展欧几里得算法。
假如我们要求一个不定方程ax+by=c的整数解(a,b,c都是非0整数)。
那么,这个不定方程有整数解。当且仅当c能被gcd(a,b)整除,
其必要性很好证明,由于以下解法,我们也能证明它的充分性。
我们假设,gcd(a,b)=g,c=C*g,那么ax+by=C*g
那么我们只要求出ax+by=g的解,就能求出原方程的解了。(在那基础上*C即可)
以下将说明
我们可以通过bx+(a%b)y=g的解,得到ax+by=g的解。这样我们就可以通过gcd(a,b)的递归过程来顺便求得不定方程的解。
假设,bx1+(a%b)y1=g①
ax2+by2=g,
已知x1,y1,求x2,y2。
方程①等价于bx1+(a-[a/b]*b)y1=g 这里[]符号是向下取整
那么bx1+(a-[a/b]*b)y1 = ax2+by2
b(x1-[a/b]*y1-y2) + a*(y1-x2) = 0
x2,y2的一个解就是{x2=y1,y2=x1-[a/b]*y1}
众所周知,正常的不定方程有无数解,
若x=X,y=Y是方程的一个解的话,
那么该不定方程的解集就是:x=X+kb,y=Y-ka,k*g为整数
于是我们得到以下代码:
int exgcd(int a,int b,int&x,int&y)
{
if(!b){x=1,y=0;return a;}
int x1,y1;
int g=exgcd(b,a%b,x1,y1);
x=y1;y=x1-a/b*y1;
return g;
}
扩展欧几里得算法,超级有用
那么重头戏来了。逆元,是个非常,非常,非常有用的东西。
几乎每道计数题,都要用到它。
在讲逆元之前,我们先来聊聊取模。
有一个东西叫同余方程,a≡b(mod p) 表示a%p=b%p。
由于中间那个符号太难打,我之后都会用第二种方法来表示啦。
我们来聊聊%的性质。
(a%p+b%p)%p = (a+b)%p
(a%p-b%p)%p = (a-b)%p
(a%p)*(b%p)%p = (a*b)%p
有一些重要的定理(摘录自别的博客):
还有一个重要的定理(费马小定理):
当p为质数时,a^(p-1) %p=1。(这个很重要,不会证也得背下来)
还有另一个较完整的定理(欧拉定理)。感兴趣的话看 http://blog.csdn.net/qq_24451605/article/details/47045279
计数问题中,有很多答案都是非常非常大的,为了避免高精度,精简答案,许多任务都要求将答案对一个数p取模。
一个最简单的问题,求n!%p,很简单吧,都会吧。
那么如果求n!/m!%p呢?(m<=n)。
好了是时候来聊一聊逆元了
逆元的定义是这样的:
假如ax % p=1,
那么最小的正整数解x就是a的逆元。(如果无解的话,那就是没有逆元了)
那么,当p不相同的话,a的逆元也不同。
所以这个p很重要。
我们来想一下,假设我们要求(s/a) %p的值,
我们设x是a的逆元(%p意义下)
那么ax %p = 1
sax % p = s(ax%p) %p = s%p
sx % p = s/a % p
所以,我们只要知道a的逆元,就可以求出(s/a) %p的值了。
当p是质数的时候:
因为a^(p-1) %p = 1(费马小定理),
所以a * a^(p-2) % p = 1
所以a^(p-2)就是a的逆元。
很完美吧。
那么有一种做法,就是当p是质数的时候,快速幂来求逆元。
那么如果p不是质数呢?
我们来想一下,ax % p =1
我们设ax = py +1(x,y均为整数),那么ax%p=1很显然吧
于是,ax-py = 1,ax + (-p) y =1
是不是很眼熟?求不定方程的整数解,扩展欧几里得算法!。
于是我们也知道,如果gcd(a,p)不为1的话,a就没有逆元了。
我们可以通过欧几里得算法,来求出该不定方程的一个解。
如果要使x为正整数且尽可能的小呢?
通过前面所说的:x=X+kb
我们将x赋为 (b+(x%b))%b 即可。(c++的取模运算会保留负号的)
那么,逆元的一部分就聊完了
我们来回顾一下
当p为质数,a^(p-2)是a的逆元,通过快速幂可求。
否则,可以用扩展欧几里得算法,设ax-py = 1,求得不定方程的解再转化,即可求逆元。
通过上面长长的赘述,
我们来聊一下组合数。
定义组合数C(n,m)为n!/m!/(n-m)!,
生活中的意义就是在n个不同的球里面无次序选m个球的方案数。
由于组合数可能非常巨大,
通常要求C(n,m)%p。
即求n!/m!/(n-m)!%p。
如果n<=1000000,m<=1000000
我们发现,除号后面的都是阶乘,
我们可以预处理阶乘,预处理阶乘的逆元,即可快速求组合数。
来思考一下,
当p是个质数的时候,求1000000个阶乘的逆元,每次快速幂很慢吧。
那么怎么快速求阶乘的逆元呢?
我们来想,假如x!的逆元是a,
那么 x! * a % p =1
于是(x-1)! * (x*a) % p =1
我们就知道了 (x-1)!的逆元 就是 x!的逆元*x%p
这样,先求出1000000!的逆元,再倒序求阶乘逆元,速度就非常快了。。
不过
有一些特殊情况:比如1000! % 300 =0
那么1000!就没有逆元了。
所以上述方法适用于p是质数且p>n,p>m的情况。
如果p不是质数??如果p<=n或p<=m??
自己去想吧,太晚了我该睡觉了。。