幼儿园数学

CAUTION
This pointless post is intended for childern under six.

To begin with,watch this video.

素数

定义

不能被 11 和它本身整除的正整数称为素数。比如,1926081719260817 就是一个素数。数学语言:

dp,d{xxN,x1,xp}d\nmid p,d\in {x\mid x\in \mathbb{N},x\neq 1,x\neq p}

在整个自然数集合中,质数的数量并不多,分布比较稀疏,对于一个足够大的整数 NN,不超过 NN 的质数大约有 NlnN\frac{N}{\ln N} 个,即每 lnN\ln N 个数中有大约 11 个质数。

Quoted frome 《算法竞赛进阶指南》

判定

根据定义,很容易想到通过 O(N)O(\sqrt{N}) 地不断试验能否被 [2,N][2,\sqrt{N}] 中的整数整除来判断 NN 是否为素数。

1
2
3
4
5
bool is_prime(int n) {
for(int i=sqrt(n);i>1;--i)
if(!n%i) return false;
return true;
}

线性筛

首先考虑一种小学数学课本上的方法:尝试 [2,x][2,x] 中的每个整数,根据素数的定义把它的倍数 2x,3x,,xNx2x,3x,\dots , x\left\lfloor\frac{N}{x}\right\rfloor 标记为合数,未被标记则为素数。在算法进行的过程中会有合数被反复标记,所以这样做的复杂度显然是非线性的。因此,我们需要让每个合数有唯一确定的表示方式:用数组 mpf 记录每个数的最小素因数,素数的 mpf 为 0,首先枚举 [2,N][2,N] 中的整数,如果未被记录最小素因数则是素数,接着去筛除它的倍数,与普通筛法不同的是,线性筛只会筛除当前枚举的数比它的最小素因数小的素数

1
2
3
4
5
6
7
8
void linearFilter(const int& n) {
for(int i=2;i<=n;++i) {
if(!mpf[i]) mpf[i]=p[++m]=i; // 未得到最小素因数,所以是素数
for(int j=1,k=n/i;j<=m;++j) // 枚举已找到的素数
if(p[j]>mpf[i]||p[j]>k) break; // 只要这个素数比当前数的最小素因数小而且乘完不过界
else mpf[i*p[j]]=p[j]; // 相乘就是合数且最小素因数是 p[j]
}
}

素因数分解

算数基本定理:任何大于 11 的正整数都能被唯一分解为若干个素数的积,即:(其中 pip_i 是递增素数)
N=i=1mpici(ciN)N=\prod\limits_{i=1}^mp_i^{c_i}(c_i\in \mathbb{N})

因此可以从小到大用 [2,N][2,\sqrt{N}] 中的整数尝试去除欲分解的数,找到一个素因数就除到不能整除为止。又根据算数基本定理,易知此法得到的一定是素因数。

1
2
3
4
5
6
7
8
void div(int n) {
for(int i=2,s=sqrt(n);i<=s;++i)
if(!(n%i)) {
p[++m]=i;
while(!(n%i)) n/=i,++c[m];
}
if(n>1) p[++m]=n,++c[m];
}

约数

定义

n,dZn,d\in \mathbb{Z},若 nmodd=0n\bmod d=0,则 dd 整除 nnddnn 的约数,nndd 的倍数,记作 dnd\mid n

算数基本定理的推论

NNN\in \mathbb{N},根据算数基本定理易知:

  • NN 的正约数个数为 i=1m(ci+1)\prod\limits_{i=1}^m(c_i+1)
  • NN 的正约数集合为 {i=1mpibi}(0bici){\prod\limits_{i=1}^mp_i^{b_i}}(0\leq b_i \leq c_i)
  • NN 所有正约数和为 i=1mj=0cipij\prod\limits_{i=1}^m{\sum\limits_{j=0}^{c_i}p_i^j}

试除法

依次尝试 [2,N)[2,\sqrt N) 中的整数,若 ddNN 的约数,则 Nd\frac{N}{d} 也是 NN 的约数。约数个数上界为 2N2\sqrt N

倍数法

试除法求 NN 个正整数的约数集合复杂度为 Θ(NN)\Theta (N\sqrt N),可以使用复杂度为 O(NlogN)O(N\log N) 的倍数法。因为 [1,N][1,N] 中有约数 dd 的整数数为 dd 的倍数,所以:

1
2
3
for(int i=1;i<=n;++i)
for(int j=1,k=n/i;j<=k;++j)
fac[i*j].push_back(i);

[1,N][1,N] 中所有整数的约数共有约 NlogNN\log N 个。

最大公约数

定义

d,m,a,bNd,m,a,b\in \mathbb{N}ddaabb 的约数,称 dda,ba,b 的公约数,其中最大的 dd 为最大公约数。记作 gcd(a,b)\gcd(a,b)mmaabb 的倍数,称 mma,ba,b 的公倍数,其中最小的 mm 为最小公倍数。记作 lcm(a,b)\text{lcm}(a,b)。此外,有几个定理:

  • a,bN gcd(2a,2b)=2gcd(a,b)\forall a,b \in \mathbb{N}\ \gcd(2a,2b)=2\gcd(a,b)
  • gcd(a,b)lcm(a,b)=ab (a,bN)\gcd(a,b)\text{lcm}(a,b)=ab\ (a,b\in \mathbb{N})
  • a{2k1kN},b{2kkN} gcd(a,b)=gcd(a,b2)a\in {2k-1\mid k\in \mathbb{N}} ,b \in {2k\mid k\in \mathbb{N}} \ \Rightarrow \gcd(a,b)=\gcd(a,\frac{b}{2})

更相减损术

a,bN abgcd(a,b)=gcd(b,ab)=gcd(a,ab)\forall a,b \in \mathbb{N}\ a\geq b \Rightarrow \gcd(a,b)=\gcd(b,a-b)=\gcd(a,a-b)

1
2
3
4
int gcd(int a,int b) {
if(a<b) a^=b^=a^=b;
return b? gcd(b,a-b):a;
}

欧几里得法

a,bN b0gcd(a,b)=gcd(b,amodb)\forall a,b \in \mathbb{N} \ b\neq 0\Rightarrow \gcd (a,b)=\gcd (b,a\bmod b)

1
2
3
int gcd(int a,int b) {
return b? gcd(b,a%b):a;
}

互素

a,b,cNa,b,c\in \mathbb{N}。若 gcd(a,b)=1\gcd(a,b)=1 则称 a,ba,b 互素。若 gcd(a,b,c)=1\gcd(a,b,c)=1 则称 a,b,ca,b,c 互素。若 gcd(a,b)=gcd(a,c)=gcd(b,c)=1\gcd(a,b)=gcd(a,c)=\gcd(b,c)=1 称为 a,b,ca,b,c 两两互素。

积性函数

定义

a,ba,b 互素,有 f(ab)=f(a)f(b)f(ab)=f(a)f(b)。那么函数 ff 为积性函数。

一些性质

  • 根据算数基本定理,f(n)=i=1mf(pici)f(n)=\prod\limits_{i=1}^{m}f(p_i^{c_i})

欧拉函数

定义

[1,N][1,N] 中与 NN 互素的数的个数称为欧拉函数,记作 φ(N)\varphi (N)

一些性质

  • 欧拉函数是积性函数,有积性函数的所有性质。
  • φ(N)=Ni=1m(11pi)\varphi (N)=N\prod\limits_{i=1}^m(1-\frac{1}{p_i})
  • n>1,[1,N]\forall n > 1,[1,N] 中与 nn 互素的数的和为 N2φ(N)\frac{N}{2}\varphi (N)
  • pnp\mid np2np^2\mid n,则 φ(n)=pφ(np)\varphi (n)=p\varphi (\frac{n}{p})
  • pnp\mid np2np^2\nmid n,则 φ(n)=(p1)φ(np)\varphi (n)=(p-1)\varphi (\frac{n}{p})
  • dnφ(d)=n\sum\limits_{d\mid n}\varphi (d)=n

求欧拉函数

在进行质因数分解的过程中就可以按照欧拉函数的公式进行计算:

1
2
3
4
5
6
7
8
9
10
int phi(int n) {
int ans=n;
for(int i=2,s=sqrt(n);i<=s;++i)
if(!(n%i)) {
ans=ans/i*(i-1);
while(!(n%i)) n/=i;
}
if(n>1) --ans;
return ans;
}

同余

定义

a,bZa,b\in \mathbb{Z}。若 amodm=bmodma\bmod m=b\bmod m,则称 a,ba,b 在模 mm 的意义下同余,记作 ab(modm).a\equiv b\pmod{m}.

{a+km0am1,a,kZ}{a+km\mid 0\leq a\leq m-1,a,k\in \mathbb{Z}} 中的数模 mm 同余,则称该集合为一个模 mm 的同余类,记作 a\overline{a}。模 mm 的同余类分别为 0,1,2,,m1\overline{0},\overline{1},\overline{2},\dots ,\overline{m-1},这 mm 个同余类构成 mm 的完全剩余系,简称完系。[1,m][1,m] 中与 mm 互素的数代表的同余类有 φ(n)\varphi (n) 个,构成 mm 的简化剩余系,简称缩系。1a,bm1\leq a,b\leq ma,ba,bmm 互素,则 abab 也与 mm 互素,再根据模运算的性质可知 abmodmab \bmod m 也与 mm 互素,即 abmodmab \bmod m 也属于 mm 的缩系,这就是缩系关于模 mm 乘法封闭。

一些定理

  • 欧拉定理:a,bN,gcd(a,b)=1aφ(n)1(modn)a,b\in \mathbb{N},\gcd (a,b)=1\Rightarrow a^{\varphi (n)}\equiv 1\pmod{n}
  • 费马小定理(即当欧拉定理中的 nn 为素数时):pp 是素数,aZ\forall a\in \mathbb{Z},有 apa(modp)a^p\equiv a\pmod{p}
  • 欧拉定理的推论:gcd(a,n)=1,bN\gcd (a,n)=1,\forall b\in \mathbb{N}ababmodφ(n)(modn)a^b\equiv a^{b\bmod \varphi (n)}\pmod{n}。根据欧拉定理的推论,计算一个乘方式对素数 pp 取模的结果时,可以先让底数对 pp 取模指数对 φ(p)\varphi (p) 取模,再计算答案。
  • Bézout 定理:a,bZ,x,y\forall a,b\in \mathbb{Z},\exists x,y,满足 ax+by=gcd(a,b)ax+by=\gcd (a,b)

Bézout 定理的证明:在欧几里得法中 b>0b>0 时有 gcd(a,b)=gcd(b,amodb)\gcd (a,b)=\gcd (b,a \bmod b),代入整理得 bx+(amodb)y=bx+(abab)y=ayb(xaby)bx+(a\bmod b)y=bx+(a-b\left\lfloor\frac{a}{b}\right\rfloor)y=ay-b(x-\left\lfloor\frac{a}{b}\right\rfloor y),所以令 x=y,y=xabyx'=y,y'=x-\left\lfloor\frac{a}{b}\right\rfloor y 就可以得到等式的解;b=0b=0 时,显然有 {x=1y=0\begin{cases}x=1\y=0\end{cases} 使等式成立。Q.E.D.Q.E.D.

扩展欧几里得

Bézout 定理的证明同时给出了 x,yx,y 所有解的计算方法,这种方法称为扩展欧几里得算法:

1
2
3
4
5
void exgcd(int a,int b,int& x,int& y) {
if(!b) {x=1,y=0;return;}
exgcd(b,a%b,x,y);
int z=x;x=y;y=z-a/b*y;
}

如果你和我一样无聊,还可以用 JavaScript 再写一遍

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
const readLine = require('readline');
const rl = readLine.createInterface(process.stdin, process.stdout);

var x, y;

rl.on('line', str => {
var [a, b] = str.split(' ').map(s => +s);
exgcd(a, b);
process.stdout.write((x - Math.floor(x / b) * b) + '');
rl.close();
});

function exgcd(a, b) {
if(!b) [x, y] = [1, 0];
else exgcd(b, a % b), [y, x] = [x - Math.floor(a / b) * y, y];
}

乘法逆元

b,mZ,gcd(b,m)=1,bab,m \in \mathbb{Z},\gcd (b,m)=1,b\mid a,则存在一个整数 xx 使得 abax(modm)\frac{a}{b}\equiv ax\pmod{m}xx 就是 bbmm 乘法逆元,记作 b1(modm)b^{-1}\pmod{m}

上式移项可得 bb11(modm)b\cdot b^{-1}\equiv 1\pmod{m}。若 pp 为素数且 b<pb<p,根据费马小定理可得 bp11(modp)bbp21(modp)b^{p-1}\equiv 1\pmod{p}\Rightarrow b\cdot b^{p-2}\equiv 1\pmod{p}。即 bp2b^{p-2}bb 的模 pp 乘法逆元。

1
2
3
4
5
6
7
8
int qpow(int a,int b,int ans=1) {
for(;b;b>>=1,a=(1LL*a*a)%P)
if(b&1) ans=(1LL*a*ans)%P;
return ans;
}
int inv(int n) {
return qpow(n,P-2);
}

对于更一般的情况,乘法逆元可以通过解一次同余方程 bx1(modm)bx\equiv 1\pmod{m} 得到。

一次同余方程

a,b,m,xZa,b,m,x\in \mathbb{Z},形如 axb(modm)ax\equiv b\pmod{m} 的方程称为一次同余方程。axb(modm)maxbax\equiv b\pmod{m} \Leftrightarrow m\mid ax-b。设 y=axbm-y=\left\lfloor\frac{ax-b}{m}\right\rfloor,于是方程可改写为 ax+my=bax+my=b。根据 Bézout 定理,上式有解当且仅当 gcd(a,m)b\gcd (a,m)\mid b。有解时,先用 exgcd 求出 ax+my=gcd(a,m)ax+my=\gcd (a,m) 的一组解,然后再乘上 bgcd(a,m)\frac{b}{\gcd (a,m)} 就得到原方程的一组解,所有模 mgcd(a,m)\frac{m}{\gcd (a,m)} 与得出的一组解同余的整数都是方程的解。

矩阵

A,B,CA,B,Cn×mn\times m 的矩阵,C=A±Bi[1,n],j[1,m],Ci,j=Ai,j±Bi,jC=A\pm B\Leftrightarrow \forall i\in [1,n],\forall j\in [1,m],C_{i,j}=A_{i,j}\pm B_{i,j}

AAn×mn\times m 矩阵,BBm×pm\times p 矩阵,则 C=ABC=ABn×pn\times p 的矩阵,第一个矩阵的列数等于第二个矩阵的行数时矩阵乘法才有意义。计算方法:

Ci,j=k=1mAi,kBk,jC_{i,j}=\sum\limits_{k=1}^{m}A_{i,k}B_{k,j}

如果没看懂请自行画图模拟。此外,矩阵乘法满足分配律和结合律,但不满足交换律。
FF1×n1\times n 矩阵,TTn×nn\times n 矩阵,F=FTF'=FT1×n1\times n 矩阵。j[1,n],Fj=k=1nFkAk,j\forall j\in [1,n],F'j=\sum\limits{k=1}^nF_kA_{k,j}。在这样的乘法运算中 FkF_kAk,jA_{k,j} 为系数累加到 FjF'_j 上。所以线性递推式一般可以写成矩阵的形式。

矩阵快速幂

因为矩阵乘法满足结合律,所以可以使用快速幂思想进行优化。比如矩阵加速地推 Fibonnacci 数列

组合数

组合数排列数的公式与性质

  • Amn=m!(mn)!A_m^n=\frac{m!}{(m-n)!}
  • Cmn=m!n!(mn)!C_m^n=\frac{m!}{n!(m-n)!}
  • Cmn=CmmnC_m^n=C_m^{m-n}
  • Cmn=Cm1n+Cm1n1C_m^n=C_{m-1}^n+C_{m-1}^{n-1}
  • i=0nCni=2n\sum\limits_{i=0}^nC_n^i=2^n

证明过程请参考《王后雄教材完全解读》

组合数计算

根据 Cmn=Cm1n+Cm1n1C_m^n=C_{m-1}^n+C_{m-1}^{n-1},组合数可以写成杨辉三角,通过 Θ(N2)\Theta (N^2) 的递推计算所有组合数。若计算某个组合数对 pp 曲目的值,可以先计算 m!m! 再计算 n!(nm)!n!(n-m)! 每项的逆元并相乘,复杂度为 O(mlogm)O(m\log m)

1
2
3
4
5
6
7
int C(int n,int m) {
if(n<m) return 0;
int ans=1;
for(int i=n-m+1;i<=n;++i) ans=(1LL*ans*i)%P;
for(int i=2;i<=m;++i) ans=(1LL*ans*inv(i))%P;
return ans;
}

本站所有原创内容采用 CC-BY-NC-ND 4.0 International 协议进行许可,附加条款亦可能应用。部分非原创内容版权为 上海アリス幻楽団黄昏フロンティア 等同人创作者所有。

笔记 算法

评论

您需要成为穿墙的邪仙以加载 Disqus。