扩展欧几里得算法学习笔记

前置芝士

exgcd 就是用来求解 $ax+by=\gcd(a,b)$ 这个方程的最小整数解的（$a,b \in \mathbb{N^*}$）。

拿到这个柿子，第一步当然是推导、化简啦。

推导过程

首先为了方便推导，令 $a\ge b$；

很明显，若 $b=0$，则 $\gcd(a,b)=a,x=1,y=0$；

对于 $b\ne 0$ 的情况：

设 $x=x^\prime,y=y^\prime$ 为 $bx+(a\%b)y=\gcd(b,a\%b)$ 的最小整数解；

$\because \gcd(a,b)=\gcd(b,a\%b)$

$\therefore ax+by=bx^\prime+(a\%b)y^\prime$

$\because a\ge b$

$\therefore a\%b=a-b\lfloor\dfrac{a}{b}\rfloor$

$\therefore ax+by=bx^\prime+(a-b\lfloor\dfrac{a}{b}\rfloor)y^\prime$

$\qquad \qquad \,\,\,= bx^\prime+ay^\prime-b\lfloor\dfrac{a}{b}\rfloor y^\prime$

$\qquad \qquad \,\,\,= ay^\prime+b(x^\prime-\lfloor\dfrac{a}{b}\rfloor y^\prime)$

而我们想要令 $x,y$ 最小，所以 $x=y^\prime,y=x^\prime-\lfloor\dfrac{a}{b}\rfloor y^\prime$。

模板代码

int exgcd(int a,int b,int &x,int &y) // x 和 y 是引用，因为我懒得写结构体…… 
{
	if(a<b) // 如果 a < b 那么交换 a 和 b 来求 
	{
		return exgcd(b,a,y,x);
	}
	if(b==0) // 如果 b 为 0，那么 gcd(a,b) = a，x = 1，y = 0 
	{
		x=1;
		y=0;
		return a;
	}
	// 其它情况 
	int tmpx,tmpy;
	int res=exgcd(b,a%b,tmpx,tmpy);
	x=tmpy;
	y=tmpx-a/b*tmpy;
	return res;
}

应用（求乘法逆元）

$exgcd$ 可以干很多事情，甚至还有 $exexgcd$……（扩展扩展欧几里得函数）

所以理解好它很重要 awa

$exgcd$ 的一个比较常见的用途是求乘法逆元。

设 $a$ 在 $\mod{b}$ 意义下的逆元是 $a^*$，那么

$a\cdot a^* \equiv 1\mod{b}$

即

$b\mid (a\cdot a^*-1)$

设 $\dfrac{a\cdot a^*-1}{b}$ 为 $y$，那么

$a^*\cdot a -1=by$

即

$a^*\cdot a-by=1$

这个方程就是 $ax+by=\gcd(a,b)$ 的形式，所以我们可以得知**$x$ 的逆元存在的必要条件是 $x$ 和模数 $p$ 互质**。

所以我们只需要使用 $exgcd$ 就可以求出逆元了。

求逆元代码如下：(P5431 【模板】乘法逆元 2)

#include <iostream>
#include <cstdio>
#include <cstring>

using namespace std;

int n;
long long p,k;
long long a[5000005],sum[2][5000005];

inline long long read()
{
	long long s=0,w=1,ch=getchar();
	while(ch<'0'||ch>'9') ch=='-'?w=-1,ch=getchar():ch=getchar();
	while(ch>='0'&&ch<='9') s=(s<<1)+(s<<3)+(ch^48),ch=getchar();
	return s*w;
}

long long exgcd(long long a,long long b,long long &x,long long &y) // x 和 y 是引用，因为我懒得写结构体…… 
{
	if(a<b) // 如果 a < b 那么交换 a 和 b 来求 
	{
		return exgcd(b,a,y,x);
	}
	if(b==0) // 如果 b 为 0，那么 gcd(a,b) = a，x = 1，y = 0 
	{
		x=1;
		y=0;
		return a;
	}
	// 其它情况 
	long long tmpx,tmpy;
	long long res=exgcd(b,a%b,tmpx,tmpy);
	x=tmpy;
	y=tmpx-a/b*tmpy;
	return res;
}

inline long long inv(long long val) // exgcd 求乘法逆元 
{
	long long x,y;
	if(exgcd(val,p,x,y)!=1)
	{
		return -1;
	}
	return (x%p+p)%p;
}

int main()
{
	scanf("%d%lld%lld",&n,&p,&k);
	for(int i=1;i<=n;i++)
	{
		a[i]=read();
	}
	sum[0][0]=1;
	sum[1][n+1]=1;
	for(int i=1;i<=n;i++)
	{
		sum[0][i]=sum[0][i-1]*a[i]%p;
	}
	for(int i=n;i>=1;i--)
	{
		sum[1][i]=sum[1][i+1]*a[i]%p;
	}
	long long summ=0,base=1;
	for(int i=1;i<=n;i++)
	{
		base=base*k%p;
		summ+=base*sum[0][i-1]%p*sum[1][i+1]%p;
		summ%=p;
	}
	printf("%lld\n",summ*inv(sum[0][n])%p);
	return 0;
}

练习题目

P2613 【模板】有理数取余

P1516 青蛙的约会