浅谈基础数论算法
数论，作为离散数学一个分支，是计算机科学的重要组成部分。在网络安全等方面有着广泛的应用。今天笔者将从基础内容出发，简要介绍最大公约数，解二元一次不定方程，解同余方程组等算法知识。
最大公约数
最大公约数，即两个整数的最大公因子。在计算机科学中，我们通常采取辗转相除的方法来计算两个数的最大公约数，代码如下。
int gcd(int a,int b)
{
if(!b)return a;
else return gcd(b,a%b);
}
这一算法的正确性源于一个简单的数论定理(gcd(a,b)=gcd(b,a-b))，由此基本定理进行简单归纳可知(gcd(a,b)=gcd(b,a\%b))。这两个定理的证明都较为基础，留给读者作为练习。
读者读到这里可能会对该算法的时间复杂度抱有一些疑惑。事实上，我们进行简单分析即可发现，在((a,b)->(b,a\%b))的过程中，我们有(a\%b leq lfloor frac{a}{2} floor)，即其中一个数字会衰减到原来的一般。由此，我们可以得知，最多经过(log_2{(a+b)})轮，两个数字其中之一就会变为(0)，此时算法也就运行结束，故该算法的时间复杂度为(O(log_2n))
二元一次不定方程
二元一次不定方程指的是形如(ax+by=c)的方程，其中(a,b,c)为参数，(x,y)为变量，取值范围均为整数。
在计算机科学中，经常会遇到求解这类方程的特解，通解，正整数解，解的数量等问题。而扩展欧几里得算法((exgcd))给出了一个求解该类方程通解的方法。
具体的算法流程是这样的，我们考虑沿用求两个数字最大公约数时的辗转相除法，迭代构造求解。
首先，在递归出口(b=0)时，(d=gcd(a,b)=a)，原方程转化为(ax=a)，此时我们取(x_0=1,y_0=0)即可得到一组特解。我们设在递归上一层时的参数为(A,B)，我们有(a=B,b=A\%B)。由刚才构造的特解可知，我们有方程(B*x_0+(A\%B)*y_0=d)成立。我们尝试构造(x,y)使得新的方程(A*x+B*y=d)。经过推导，我们可以得到以下表达式：
代入验证不难发现这组新的((x,y))一定可以使得(A*x+B*y=d)成立。因此，我们只需要回溯即可递推出原方程的一组特解，设其为(X_0,Y_0)，我们可以证明，对于原方程(ax+by=1)，通解形式如下：
证明较为复杂，读者可以自行查看相关资料。
如下程序展示了求解二元一次不定方程非负整数解的方法，供以参考。
int X,Y;
int exgcd(int a,int b)
{
if(!b){X=1;Y=0;return a;}
int d=exgcd(b,a%b),t=X;X=Y;Y=t-(a/b)*Y;
return d;
}
int main()
{
int a=read(),b=read(),c=read();//ax+by=c
int d=exgcd(a,b);//d is the gcd of a and b
if(c%d!=0)printf("No solution");
else
{
X*=c/d;Y*=c/d;
int p1=b/d,p2=a/d;//x=k*X+t*p1 y=k*Y-t*p2
X=(X%p1+p1)%p1;Y=(c-a*X)/b;
if(Y<0)printf("There is no positive integer solution");
else printf("%d %d
",X,Y);
}
return 0;
}
扩展欧几里得算法还可以同时计算正整数解数量，(x)最大的正整数解，(y)最大的正整数解等问题，限于篇幅限制，本文不再赘述，具体可以参考这道例题https://www.luogu.com.cn/problem/P5656。
利用这一算法，我们可以解决许多其他问题。比如求整数(a)在模(p)意义下的除法逆元(x)，即(ax equiv 1 pmod{p})，我们只需要把该式转化为(ax+kp=1)即可。此外，该算法在接下来介绍的解同余方程组中也有应用。
同余方程组
同余方程组指的是对于许多条形如(xequiv a pmod{p})的方程，去计算它们解集的交集。
我们考虑对于两条方程
等价于：
联立消元得：
这显然是一个关于(k_1,k_2)的二元一次不定方程，套用扩展欧几里得算法即可得到一组(k_1,k_2)的特解。
然后根据中国剩余定理，我们可以得知，(x equiv a_1+k_1*p_1 pmod{lcm(p_1,p_2)})，这样我们就将两个方程合并成了一个方程。只需要不断沿用此方法，即可将方程组两两合并成一个方程，得到方程组的通解。
代码如下：
ll ksc(ll x,ll y,ll p){return ((x*y-(ll)(((ldb)x/p*y))*p)%p+p)%p;}//防止溢出的乘法
ll X,Y;
ll exgcd(ll a,ll b)
{
if(!b){X=1,Y=0;return a;}
ll d=exgcd(b,a%b),t=X;X=Y;Y=t-(a/b)*Y;
return d;
}
ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}
ll lcm(ll a,ll b){return (a/gcd(a,b))*b;}
int main()
{
// x=a (mod p)
ll n=read(),p,a;
for(ll i=1;i<=n;i++)
{
ll tp=read(),ta=read();
if(i==1){p=tp;a=ta;continue;}
ll np=lcm(p,tp),d=exgcd(p,tp);
X=ksc(X,(ta-a)/d,tp/d);
a=(ksc(X,p,np)+a)%np;p=np;//合并结果
}
printf("%lld",a);
return 0;
}