同余系基本定理1

· · 个人记录

我与同余系的第一次相逢，是在一个寒冷的冬日。

胡扯完毕，大家放轻松。

\tiny \text{(2020.4.11) 两年天坑终于更完了!}

-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-

1.同余系基本运算

首先讲一下同余系的概念,对于第一次接触除了经典实数体系之外的同学,可能有点难理解。

大家习惯的是实数域 ? ,如果您告诉我您习惯 ? 您可以马上跳过这一节。

这个年代,有很多很多的计数题目,这类题目的答案往往是天文数字。

在遥远的过去,毒瘤出题人往往会要求选手写高精度来求出确切的答案,于是计数题一直没有普及。

直到同余系普及,计数题目便不再依赖大数计算,得以普及并迅速发展。于是就出现了新的毒瘤。

或许你曾见过计数题带有模数,如 P4017 最大食物链计数 , P5520 [yLOI2019] 青原樱 , P5888 传球游戏 等等。

做这种题的时候,一种容易想到的思路是 : 先求出确切答案,然后再取模。

当然这是相当搞笑的行为,如果你第一次见到对答案取模,去询问学长或者老师,他们多半会告诉你:

(a+b)\bmod p=(a\bmod p)+(b\bmod p)\bmod p (a*b)\bmod p=((a\bmod p)*(b\bmod p))\bmod p

太过通俗,不证。

也就是说,涉及加法和乘法时,取模的顺序不会影响答案。

这启发我们建立同余系(模p意义下):

\color{blue}\text{定义:}

这样,同余系中的数很小,就不需要涉及大数运算了。

当然同余系存在的理由并不止于计数题,它本身就是一个优美而深奥的世界,可以导出许多有用的结论。

那么,就让我们开始吧!

部分计数题目,只涉及加法和乘法,我们直接使用上面介绍的同余系就可以解决。

可能你会好奇 : 除法去哪了?这里都是整数,除出分数怎么办?

我们当时是怎么定义除法的呢?

根据除法是乘法的逆运算,我们要找到一个数,使得其能“抵消”乘法的影响,这就是倒数。

就是构造出\frac{1}{a},使得a*\frac{1}{a}=1,然后乘以\frac{1}{a}这个操作就等价于除以a.

类似地,我们可以\color{blue}\text{定义}同余系内的乘法逆元(倒数)a^{-1}是 : 满足ab=1\pmod p的b

根据这个同余系的定义,这个a^{-1}是整数。没错,不是分数是整数。

这就是同余系的魅力了,无论在经典数系里面怎样的复杂,只要能在同余系中表示,就一定是整数。

$a=0$时$a^{-1}$不存在,正如经典数系中不能除以$0$. 其他情况下,$a^{-1}$如果存在,则是唯一的。 $\color{blue}\text{证明:}$ 假设有$ax=1,ay=1

则有:xay=(xa)y=y,同时有xay=x(ay)=x,所以有x=y.

这个逆元怎么求呢?可以暴力枚举1...(p-1)来尝试,复杂度O(p).

显然这是很糟糕的复杂度,至于如何快速求解,后面再说。

现在四则运算都齐了,而且是封闭的,你可以像小学数学一样随意使用同余系了。

许多在实数系里面成立的东西在同余系里也成立,这里不再赘述了,请读者自行探究。

这里要求模数p是素数,由于素数有很多良好的性质,一般题目取模都是用素数,这个定理应用也最广泛。

那么,a^{-1}=a^{p-2},写个快速幂就能够求逆元了。

求逆元快速幂合一代码:

ll powM(ll a,int t=mod-2)
{
  ll ret=1;
  while(t){
    if (t&1)ret=ret*a%mod;
    a=a*a%mod;t>>=1;
  }return ret;
}

遗憾的是,逆元并非像实数系中那样总是存在,来看一下其存在的条件。

这玩意和逆元有什么关系呢?

注意到,当\gcd(a,b)=1的时候,方程变为ax+by=1

ax+by=1$等式两边模b得到 : $ax=1\pmod p

这正是逆元的定义式!

而且根据斐蜀定理，gcd(a,b)>1时这个方程无解。

否则gcd(a,b)=1,此时求出的x就是a^{-1}

斐蜀定理还可以扩展到更多元的情况 : P4549 【模板】裴蜀定理

P1082 同余方程 (注意模数不一定是质数)

EXgcd,即扩展欧几里得算法。

可以求出 ax+by=\gcd(a,b) 的一组解。(其中a,b已知)

这个东西的实现类似于数学中的归纳法，对初学者来说有些难理解。

\large{ax+by=\gcd(a,b)}

我们设a_2=b;\ b_2=a\%b

根据普通欧几里得算法可以得到gcd(a,b)=gcd(a_2,b_2)

假设我们知道 \large{a_2x+b_2y=gcd(a,b)} 的解 \large{x_2,y_2} 。

把a_2=b;\ b_2=a\%b代入回去。

得到bx_2+(a\%b)y_2=\gcd(a,b)

bx_2+(a-b\left\lfloor{a/b}\right\rfloor)y_2=\gcd(a,b) bx_2+ay_2-b\left\lfloor{a/b}\right\rfloor{y_2}=\gcd(a,b) 比对 $ax+by=\gcd(a,b)

得到 x=y_2;y=x_2-\left\lfloor{a/b}\right\rfloor{y_2}

也就是说我们知道a_2x+b_2y=\gcd(b,a\%b)的解之后可以O(1)推导出ax+by=\gcd(a,b)的解。

这是个嵌套的形式，对于a_2x+b_2y=\gcd(b,a\%b)，我们也转化一下，再转化，还转化……

但是这并不会无限的进行下去。b==0时,根据普通欧几里得算法，得到gcd(a,0)=a

得 ax=a，直接令x=1;y=0即为一组合法解。

这样的话，就能解上一个方程，上上个，上上上个也迎刃而解。

很明显，EXgcd的复杂度与普通的欧几里得算法相同，均为O(\log a)

我们来手玩一下加深印象。

求出 {7x+3y=gcd(3,7)} 的一组解

最终得到 x=1;y=-2

带入{7x+3y=gcd(3,7)}得到{7*1+3*(-2)=1},一切正常。

代码实现: 请仔细查看引用关系。

void exgcd(ll a,ll b,ll &x,ll &y)
{
  if (a==1&&b==0)
    {x=1;y=0;return ;}
  exgcd(b,a%b,y,x);
  y-=x*(a/b);
}

那么问题来了，这个东西和逆元有啥关系？

根据上文的推理,当\gcd(a,p)=1时,我们得到的x即为所求的逆元。

Exgcd并不基于同余系,所以可能求出来负数,模一下变成正数就好。

P5656 【模板】二元一次不定方程(exgcd)

这个毒瘤模板还要求给出最大最小解以及解的个数……可以加深对这个算法的理解。

我们的扩展欧几里得只能够处理ax+by=\gcd(a,b)的情况,而现在要求解ax+by=c的一般情况。

根据斐蜀定理,必须要满足c是d=\gcd(a,b)的倍数,否则无解。

我们令a'=a/d,b'=b/d,则有\gcd(a',b')=1,这就是我们的经典逆元方程a'x'+b'y'=1了。

直接扩欧之后,特解x_*,y_*之后,通解则是\begin{cases}x_t=x_*+tb'\\y_t=y_*-ta'\end{cases}

充分性是容易验证的,至于必要性,不证。

将一组 x,y 分别乘上 c/d 即可得到原方程的一组解ax_*+by_*=c。

其最小正整数解是容易求的,直接取模意义下的解即可。根据相应的等式可求得另一个元。

容易发现x变小时y增大,则x取最小正整数解时y是最大解,反之亦然。

可以据此判定有没有正整数解。注意可能模出0,这也是题意不自然的地方。

解的数量就是\dfrac{x_{\max}-x_{\min}}{b}+1。

#include<cstdio>
#define ll long long
#define pf printf
using namespace std;
inline int read(){
  int X=0;char ch=0;
  while(ch<48||ch>57)ch=getchar();
  while(ch>=48&&ch<=57)X=X*10+(ch^48),ch=getchar();
  return X;
}
int gcd(int a,int b)
{return !b ? a : gcd(b,a%b);}
void exgcd(ll a,ll b,ll &x,ll &y)
{
  if (b==0){x=1;y=0;return ;}
  exgcd(b,a%b,y,x);y-=x*(a/b);
}
int T;
ll a,b,c,d,x,y;
int main()
{
  int T=read();
  for (int i=1;i<=T;i++){
    a=read();b=read();c=read();
    d=gcd(a,b);a/=d;b/=d;
    if (c%d){puts("-1");continue;}
    exgcd(a,b,x,y);
    x*=c/d;y*=c/d;
    ll xl,xr,yl,yr;
    xl=(x%b+b)%b;if (!xl)xl=b;
    yr=(c-a*d*xl)/(b*d);
    yl=(y%a+a)%a;if (!yl)yl=a;
    xr=(c-b*d*yl)/(a*d);
    if (yr<=0){
      pf("%lld %lld\n",xl,yl);
    }else 
      pf("%lld %lld %lld %lld %lld\n",(xr-xl)/b+1,xl,yl,xr,yr);
  }return 0;
}

首先,这是对答案取模,而非对数列取模,看错题可能导致您神游到三里屯去。

众所周知,斐波那契数列的增长速度是指数级的,我们只需要考虑前O(\log k)项即可。

设F[n]为斐波那契数列。

然后根据递推(转移矩阵)的结合律,能得到G[n]=a*F[n-1]+b*F[n-2].

然后针对每一位,能得到形如F[n-1]x+F[n-2]y=k的不定方程。

而k不可能在数列中出现两次,每个位置的解没有交集,所以直接讲每个位置的方程的解的个数加起来即可。

似乎还是没有用到脑子,怎么办啊……

#include<cstdio>
#define ll long long
#define pf printf
using namespace std;
int gcd(int a,int b)
{return !b ? a : gcd(b,a%b);}
void exgcd(ll a,ll b,ll &x,ll &y)
{
  if (b==0){x=1;y=0;return ;}
  exgcd(b,a%b,y,x);y-=x*(a/b);
}
int T;
ll calc(ll a,ll b,ll c)
{
  ll d=gcd(a,b),x,y;a/=d;b/=d;
  if (c%d)return 0;
  exgcd(a,b,x,y);
  x*=c/d;y*=c/d;
  ll xl,xr,yl,yr;
  xl=(x%b+b)%b;if (!xl)xl=b;
  yr=(c-a*d*xl)/(b*d);
  if (yr<0)return 0;
  yl=(y%a+a)%a;if (!yl)yl=a;
  xr=(c-b*d*yl)/(a*d);
  return (xr-xl)/b+1;
}
ll k,F[105];
int main()
{
  scanf("%lld",&k);
  F[0]=F[1]=1;
  ll ans=0;
  for (int i=2;F[i-1]+F[i-2]<=k;i++){
    F[i]=F[i-1]+F[i-2];
    ans=(ans+calc(F[i-1],F[i-2],k))%1000000007;
  }printf("%lld\n",ans);
  return 0;
}

这部分需要一定的附加知识,看不懂建议跳过。

有的时候模数不一定是质数,这时候就要用到欧拉定理。这是费马小定理的扩展。

x^{\varphi(m)}=1\pmod m

其中\varphi(m)=\sum\limits_{i=1}^m[i\perp m],是欧拉函数。

证明和费马小定理类似,考虑构造与 m 互质的数的集合 S,易得|S|=\varphi(m).

将每个数乘上x得到新集合 S',由于 x\perp m,所得 S' 内的元素每个仍然与m互质,所以和原来的集合S相同。

考虑两个集合的乘积,设 S 内元素的乘积为 c,则 S' 内元素乘积为 cx^{\varphi(m)},则有c=cx^{\varphi(m)}\pmod m

由于c\perp m,方程两边可以消去c,就得到x^{\varphi(m)}=1\pmod m

令m=p可得\varphi(p)=p-1,即x^{p-1}=1\pmod p,正是费马小定理。

可以用于求逆元 : x^{\varphi(m)-1}=x^{-1}\pmod m,但是真正的用处是降幂,这会在后续文章中讲到。

P3811 【模板】乘法逆元

我们要预处理出[1,n]的逆元\pmod p(质数),直接使用上述的某一种方法大力求，复杂度为O(nlogn)。

下面介绍一些复杂度为O(n)的算法。

int inv[MaxN];
void Init()
{
  inv[1]=1;
  for (int i=2;i<=n;i++)
    inv[i]=1ll*(mod-mod/i)*inv[mod%i]%mod;
}
ll fac[MaxN],ifac[MaxN];
ll C(int n,int m)
{return fac[n]*ifac[m]%mod*ifac[n-m]%mod;}
void Init()
{
  fac[0]=1;
  for (int i=1;i<=n;i++)
    fac[i]=fac[i-1]*i%mod;
  ifac[n]=inv(fac[n]);
  for (int i=n;i;i--)
    ifac[i-1]=ifac[i]*i%mod;
}

一道比较模板的题目 : P1641 [SCOI2010]生成字符串

P5431 【模板】乘法逆元2

保证模数为质数。

观察造阶乘法的核心思想 : n^{-1}=\dfrac{(n-1)!}{n!}=fac[n-1]*ifac[n]

实际上是前缀逆乘上前缀积罢了。

考虑维护前缀积s[n]和前缀逆inv[n],方法同上。

然后使用(a[n])^{-1}=inv[n]*s[n-1]即可。

复杂度是O(n+\log p)的,在某些时候能派上用场。

注意实际应用中,如果有0出现需要特判。

#include<cstdio>
#define ll long long
#define MaxN 5000500
using namespace std;
inline int read(){
  register int X=0;
  register char ch=0;
  while(ch<48||ch>57)ch=getchar();
  while(ch>=48&&ch<=57)X=X*10+(ch^48),ch=getchar();
  return X;
}
int mod;
ll powM(ll a,int t=mod-2)
{
  ll ret=1;
  while(t){
    if (t&1)ret=ret*a%mod;
    a=a*a%mod;t>>=1;
  }return ret;
}
ll inv[MaxN],s[MaxN],k;
int n,a[MaxN];
int main()
{
  scanf("%d%d%lld",&n,&mod,&k);
  s[0]=1;
  for (int i=1;i<=n;i++)
    s[i]=s[i-1]*(a[i]=read())%mod;
  inv[n]=powM(s[n]);
  for (int i=n;i;i--)
    inv[i-1]=inv[i]*a[i]%mod;
  ll ans=0,buf=1;
  for (int i=1;i<=n;i++){
    buf=buf*k%mod;
    ans=(ans+buf*inv[i]%mod*s[i-1])%mod;
  }printf("%lld",ans);
  return 0;
}

欲知后事如何,请见 同余系基本定理2