数论学习笔记

素数:

#include<iostream>#include<cstdio>#include<cstring>using namespace std;int p[10000000],a[100000001];int main(){    int i,t=0,j,n;    scanf("%d",&n);    a[1]=1;    for (i=2;i<=n;++i)      {        if (!a[i])          p[++t]=i;        for (j=1;j<=t&&i*p[j]<=n&&(a[i*p[j]]==0||a[i*p[j]]>p[j]);++j)          {            a[i*p[j]]=p[j];            if (!i%p[j]) break;          }      }}

快速幂:

int pow(int x,int a,int p){    int t;    if (a==0) return 1%p;    if (a==1) return x%p;    t=pow(x,a/2,p);    if (a%2) return ((t*t)%p*x)%p;    else return (t*t)%p;}

欧拉函数:

欧拉函数ϕ(n)是指1~n-1内与n互质的数的个数，通常用公式

ϕ(n)=n∗(1−1p1)∗(1−1p2)∗...∗(1−1pk)

m=floor(sqrt(n));phi[n]=n;for (int i=2;i<=m;++i)  if (n%i==0)    {      phi[n]=phi[n]/i*(i-1);      while (n%i==0) n/=i;    }if (n>1) phi[n]=phi[n]/n*(n-1);int phi[maxn];void phin(int n){  for (int i=2;i<=n;++i) phi[i]=0;    phi[1]=1;  for (int i=2;i<=n;++i) if (!phi[i])    for (int j=i;j<=n;j+=i) {      if (phi[j]) phi[j]=j;        phi[j]=phi[j]/i*(i-1);    }}

欧拉定理:
设gcd(x,p)=1 (若x<p，则称x为p的剩余系)则有

xϕ(p)=1(modp)

特别地 当p为质数时 ϕ(p)=p−1 则有

xp−1=1mod(p)

最小公倍数类：

gcd(a,b)∗lcm(a,b)=a∗b

欧几里德算法(运算次数最多为fib[n]与fib[n-1],辗转相除):

int gcd(int a,int b){   if (a%b==0) return b;   else return gcd(b,a%b);}

高精度(Stein 算法):
gcd(x,y)=gcd(x,y/2)(x−>odd,y−>even)
gcd(x,y)=2∗gcd(x/2,y/2)(x,y−>even)
gcd(x,y)=gcd(x,y−x)(x,y−>odd)

拓展欧几里德算法，用于求ax+by=gcd(a,b))的最小解(|x|+|y|最小):

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

同余类:

线性同余方程:
ax=b(modp)=>ax−py=b=>若gcd(a,p)|b 则有扩展欧几里德算法x=x∗b/(gcd(a,p))
特别地当b=1且gcd(a,p)=1时，a与x互为模p意义下的逆元
求逆元:
基于扩展欧几里德：
ax−py=1 求得x，注意此时x并非最小正整数解x=(x+p)%p

int inv(int a,int p){    int d,x,y;    exgcd(a,p,d,x,y);    if (d==1) return (x+p)%p;    else return -1;}

基于欧拉定理
aϕ(p)=1(modp)=>aϕ(p)−1∗a=1(modp)=>x=aϕ(p)−1(modp)

线性同余方程组：
中国剩余定理(模数互质):
x=b[1](modm[1])
x=b[2](modm[2])
…
x=b[3](modm[3])
设M=Πni=1m[i]则有

x=∑i=1n(exgcd(M/m[i]%m[i],m[i])+m[i])%m[i]∗b[i]∗M/m[i]%M

#include<cstdio>using namespace std;int b[100],m[100];int n,x,y,N;void exgcd(int a,int b){    int t;    if (a%b==0)      {        x=0;        y=1;      }    else      {        exgcd(b,a%b);        t=x;        x=y;        y=t-(a/b)*y;      }}int main(){    int i,now,ans;    scanf("%d",&n); N=1; ans=0;    for (i=1;i<=n;++i)      {        scanf("%d%d",&b[i],&m[i]);        N=N*m[i];      }    for (i=1;i<=n;++i)      {        now=N/m[i];        now=now%m[i];        exgcd(now,m[i]);         x=(x%m[i]+m[i])%m[i];        x=(x*(N/m[i])*b[i])%N;         ans=(ans+x)%N;        }    printf("%d\n",ans);}

模数不互质:考虑合并方程x=b1(modm1)x=b2(modm2)
用扩展欧几里德求出x∗m1+b1=y∗m2+b2的解
使b=x∗m1+b1,m=lcm(m1,m2)得到方程x=b(modm)

a^x=b (mod p):

BSGS法(p为质数或p,a互质,此时下式中p−1为ϕ(p)):
设m=(√p),记am关于p的逆元为v=ap−1−m
暴力枚举 1≤x≤m 分别模p的值，记作e，此时考虑m+1≤x′≤2m−1
若有ax′=b(modp)则有ei∗am=b(modp)=>ei=b∗a−m(modp)=>ei=b∗v(modp)
对于im≤x≤im+m−1循环i，迭代b即可。

int solve(int a,int b,int p){    int m,v,e=1,i;    m=(int)(sqrt(p));    v=inv(pow(a,m,p),p); /*基于扩展欧几里德*/    v=pow(a,n-1-m,p);  /*基于费马小定理*/     map <int,int> x; x[1]=0;    for (i=1;i<=m;++i)      {        e=(e*a)%p;        if (!x.count(e)) x[e]=i;       }    for (i=0;i<m;++i)      {        if (x.count(b)) return i*m+x[b];        b=(b*v)%n;      }    return -1;}

扩展BSGS法(a,p不互质):
设g=gcd(gcd(a,b),p),原方程写为(a′g)x=b′g(modp′g)=>(a′)x∗gx=b′g(modp′g)=>约去g得(a′)x∗gx−1=b′(modp′),则用BSGS法即可

组合数及组合数取模类:

Cmnn!/m!(n−m)!
Cmn=Cm−1n−1+Cmn−1
Cmn=(n/m)∗Cm−1n−1=Cm−1n∗((n−m+1)/m)

组合数取模
当n,m较小p为质数:直接上逆元
当n,m较大p为质数:Lucas定理
Lucas定理:
n=nk∗pk+nk−1∗pk−1+...+n1∗p+n0
m=mk∗pk+mk−1∗pk−1+...+m1∗p+m0 则:

Cmn=Πki=0Cmini(modp)

再上逆元 (mi!(ni−mi)!)−1=(mi!(ni−mi)!)p−2

m次方求和公式：

设

Sm(n)=∑k=1nkm=1m+2m+3m+...+nm

普通求法:

nm+1−(n−1)m+1=anm+bnm−1+cnm−2+...+n

(n−1)m+1−(n−2)m+1=a(n−1)m+b(n−1)m−1+c(n−1)m−2+...+n−1

(n−2)m+1−(n−3)m+1=a(n−2)m+b(n−2)m−1+c(n−3)m−2+...+n−2

…

1m+1−0m+1=a1m+b1m−1+c1m−2+...+1

方程左右同时累加,再移项:

Sm(n)=nm+1−b∗Sm−1(n)−c∗Sm−2(n)−...−S1(n)a

伯努利数法:

Sm(n)=1m+1∑k=0mCkm+1Bknm+1−k

伯努利数：易得出Sm(0)=0易得出递推公式：

B0=1

Bn=−1n+1∑k=0n−1Ckn+1Bk

原根类:

定义:g1,g2,...,gϕ(p)能遍历p的剩余系,则称g为p的原根
原根数目为ϕ(ϕ(n))
原根性质 若ga=gb(modp)=>a=b(modϕ(p))

如何求原根，枚举判断m，是否有d|phi(p)且md=1(modp)则m不是原根
证明，若m是原根，遍历序列中一定以1结尾(欧拉定理)，若其中存在另外一个1，那么原循环节一定是新循环节的整数倍，则d|ϕ(p),即只验证ϕ(p)约数即可。

利用原根求xa=b(modp)的解:
设p的原根为g,x=gy,b=gz,利用BSGS求出z
则有gay=gz(modp)=>ay=z(modϕ(p)),利用扩展欧几里德算法，求出y
在用快速幂求得x即可

高斯消元：

#include<iostream>#include<cstring>#include<cstdio>#include<algorithm>using namespace std;double a[1000][1000],b[1000],x[1000];int n;int main(){    int i,j,k;    double xi,bz;    scanf("%d",&n);    for (i=1;i<=n;++i)       {        for (j=1;j<=n;++j)          {            scanf("%lf",&xi);            a[i][j]=xi;          }        scanf("%lf",&b[i]);      }    for (k=1;k<=n;++k)      {        if (a[k][k]==0)          {            for (j=k+1;j<=n;++j)              if (a[j][k]!=0)                {                  for (i=1;i<=n;++i)                    swap(a[j][i],a[k][i]);                  swap(b[j],b[k]);                  break;                }          }        for (j=k+1;j<=n;++j)          {            bz=a[j][k]/a[k][k];            for (i=1;i<=n;++i)              a[j][i]-=a[k][i]*bz;            b[j]-=b[k]*bz;          }      }    for (k=n;k>=1;--k)      {        x[k]=b[k]/a[k][k];        for (j=k-1;j>=1;--j)          {            b[j]-=x[k]*a[j][k];            a[j][k]=0;          }      }    double a=5-(5/0.0*0.0);    cout<<a<<endl;    for (i=1;i<=n-1;++i)      printf("%.3f ",x[i]);    printf("%.3f\n",x[n]);}