做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？？
自己去想吧，太晚了我該睡覺了。。