Acwing《算法基础课》第4章 数学知识

    科技2024-04-22  101

    Acwing《算法基础课》第4章 数学知识

    文章目录

    Acwing《算法基础课》第4章 数学知识质数判断质数分解质因数筛质数朴素筛法埃氏筛法线性筛法 约数求所有约数约数个数定理约数之和定理例子: 12 = 2 2 × 3 1 12=2^2\times3^1 12=22×31 12 12 12的约数有 1 , 2 , 3 , 4 , 6 , 12 1, 2, 3, 4, 6, 12 1,2,3,4,6,12,约数之和为 28 28 28,根据公式计算同样是 ( 2 0 + 2 1 + 2 2 ) × ( 3 0 + 3 1 ) = 28 (2^0+2^1+2^2)\times(3^0+3^1)=28 (20+21+22)×(30+31)=28个 最大公约数 欧拉函数求单个欧拉函数筛法求多个欧拉函数 快速幂扩展欧几里得算法中国剩余定理高斯消元组合数递归法预处理逆元Lucas定理分解质因数法 卡特兰数容斥原理NIM游戏取石子取石子(升级版)


    质数

    判断质数

    模板:

    bool is_prime(int x) { if (x < 2) return false; for (int i = 2; i <= x / i; i ++ ) // 注意循环条件 if (x % i == 0) return false; return true; }

    说明:

    试除法判断条件的选择 i < n和i <= n / 2时间复杂度都是 O ( n ) O(n) O(n),过高i * i <= n虽然时间复杂度是 O ( n 1 2 ) O(n^{\frac{1}{2}}) O(n21),但i * i可能会溢出因此最好的判别条件是i <= n / i,时间复杂度是 O ( n 1 2 ) O(n^{\frac{1}{2}}) O(n21)

    分解质因数

    模板:

    void divide(int x) { for (int i = 2; i <= x / i; i ++ ) // 注意循环条件 if (x % i == 0) { int cnt = 0; while (x % i == 0) x /= i, cnt ++ ; cout << i << ' ' << cnt << endl; } if (x > 1) cout << x << ' ' << 1 << endl; // 至多存在一个大于sqrt(n)的质因子 }

    说明:

    试除法数学知识: n n n最多包含一个大于 n \sqrt{n} n 的质因子。例如 6 = 2 × 3 6=2\times 3 6=2×3 6 = 2.44949 \sqrt{6}=2.44949 6 =2.44949,存在一个 > 6 >\sqrt{6} >6 的质因子 3 3 3判别质数用i <= n / i条件最好时间复杂度 O ( log n ) O(\text{log}n) O(logn),最坏时间复杂度 O ( n 1 2 ) O(n^{\frac{1}{2}}) O(n21)

    筛质数

    朴素筛法

    模板:

    int primes[N], cnt; // primes[]存储所有素数(静态链表) bool st[N]; // st[x]存储x是否被筛掉 void get_primes(){ for(int i = 2; i <= n; i++){ if (st[i]) continue; if(is_prime[i]) primes[cnt++]=i; //把素数存起来 for(int j = i; j <= n; j += i) //不管是合数还是质数,都用来筛掉后面它的倍数 st[j]=true; } }

    说明:

    思想:素数的倍数一定不是素数缺点:没必要用合数剔除时间复杂度是 O ( n log n ) O(n\text{log}n) O(nlogn)

    埃氏筛法

    模板:

    int primes[N], cnt; // primes[]存储所有素数 bool st[N]; // st[x]存储x是否被筛掉 void get_primes(int n) { for (int i = 2; i <= n; i ++ ) { if (!st[i]) { primes[cnt++] = i; // i是当前可用的最小质数,保存到primes中 for (int j = i + i; j <= n; j += i) st[j] = true; // 素数的倍数一定不是素数 } } }

    说明:

    时间复杂度是 O ( n loglog n ) O(n\text{loglog}n) O(nloglogn)核心思想:只用质数剔除特点:每次发现的第1个非标记数,一定是质数

    线性筛法

    模板:

    int primes[N], cnt; // primes[]存储所有素数 bool st[N]; // st[x]存储x是否被筛掉 void get_primes(int n) { for (int i = 2; i <= n; i ++ ) { if (!st[i]) primes[cnt ++ ] = i; for (int j = 0; primes[j] <= n / i; j ++ ) { st[primes[j] * i] = true; // 用最小质因子去筛合数primes[j] * i if (i % primes[j] == 0) break; // 若prime[j]是i的最小质因子,则prime[j+1] * i的最小质因子依旧是prime[j] } } }

    说明:

    核心思想:每个合数必有一个最小素因子,用这个因子筛掉合数理解 当i % primes[j] != 0时,说明此时遍历到的primes[j]不是i的质因子,那么只可能是此时的primes[j] < i的最小质因子,所以primes[j] * i的最小质因子就是primes[j]当有i % primes[j] == 0时,说明i的最小质因子是primes[j],因此primes[j] * i的最小质因子也就应该是prime[j],之后接着用st[primes[j+1] * i] = true去筛合数时,就不是用最小质因子去更新了,因为i有最小质因子primes[j] < primes[j+1],此时的primes[j+1]不是primes[j+1] * i的最小质因子,此时就应该退出循环,避免之后重复进行筛选 时间复杂度是 O ( n ) O(n) O(n)

    约数

    求所有约数

    模板:

    vector<int> get_divisors(int x) { vector<int> res; for (int i = 1; i <= x / i; i ++ ) if (x % i == 0) { res.push_back(i); if (i != x / i) res.push_back(x / i); } sort(res.begin(), res.end()); return res; }

    说明:

    试除法求所有约数注意区别因数分解遍历范围是1~sqrt(x),每次找到约数时,把自身和另一个同时放入最后需要排序,但实际上消耗的时间不多,int最大值才有大约1500个约数

    约数个数定理

    n = p 1 c 1 × p 2 c 2 ∗ × . . . × p k c k n = p_1^{c_1}\times p_2^{c_2} *\times ... \times p_k^{c_k} n=p1c1×p2c2×...×pkck的约数个数等于 ( c 1 + 1 ) ( c 2 + 1 ) . . . ( c k + 1 ) (c_1+1)(c_2+1)...(c_k+1) (c1+1)(c2+1)...(ck+1) 模板:

    unordered_map<int, int> primes; // 用哈希表保存质数的指数 // 质数分解 for (int i = 2; i <= x / i; i++) while(x % i == 0) { x /= i; primes[i]++; } if (x > 1) primes[x]++; // 约数个数定理 LL res = 1; for (auto elem : primes) res = res * (elem.second + 1);

    说明:

    数学原理: 一个整数能唯一展开成若干个质数乘积的形式 n = p 1 c 1 × p 2 c 2 × . . . × p k c k n = p_1^{c_1}\times p_2^{c_2} \times ... \times p_k^{c_k} n=p1c1×p2c2×...×pkck整数相乘等价于对应质指数相加 n 1 × n 2 = ( p 1 a 1 × p 2 a 2 × . . . × p k a k ) × ( p 1 b 1 × p 2 b 2 × . . . × p k b k ) = p 1 a 1 + b 1 × p 2 a 2 + b 2 × . . . × p k a k + b k n_1 \times n_2= (p_1^{a_1}\times p_2^{a_2} \times ... \times p_k^{a_k}) \times (p_1^{b_1}\times p_2^{b_2} \times ... \times p_k^{b_k})= p_1^{a_1+b_1}\times p_2^{a_2+b_2} \times ... \times p_k^{a_k+b_k} n1×n2=(p1a1×p2a2×...×pkak)×(p1b1×p2b2×...×pkbk)=p1a1+b1×p2a2+b2×...×pkak+bk 例子: 12 = 2 2 × 3 1 12=2^2\times3^1 12=22×31 12 12 12的约数有 1 , 2 , 3 , 4 , 6 , 12 1, 2, 3, 4, 6, 12 1,2,3,4,6,12 6 6 6个,根据公式计算同样是 ( 2 + 1 ) × ( 1 + 1 ) = 6 (2+1)\times(1+1)=6 (2+1)×(1+1)=6

    约数之和定理

    s u m = ∏ i = 1 k ∑ j = 0 c i p i j = ( p 1 0 + p 1 1 + . . . + p 1 c 1 ) × . . . × ( p k 0 + p k 1 + . . . + p k c k ) sum=\prod_{i=1}^k{\sum_{j=0}^{c_i}{p_{i}^{j}}}=(p_1^0 + p_1^1 + ... + p_1^{c_1})\times...\times(p_k^0 + p_k^1 + ... + p_k^{c_k}) sum=i=1kj=0cipij=(p10+p11+...+p1c1)×...×(pk0+pk1+...+pkck)

    模板:

    LL res = 1; for (auto elem : primes) { int p = elem.first, a = elem.second; LL sum = 1; while(a--) sum = sum * p + 1; res *= sum; }

    说明:

    例子: 12 = 2 2 × 3 1 12=2^2\times3^1 12=22×31 12 12 12的约数有 1 , 2 , 3 , 4 , 6 , 12 1, 2, 3, 4, 6, 12 1,2,3,4,6,12,约数之和为 28 28 28,根据公式计算同样是 ( 2 0 + 2 1 + 2 2 ) × ( 3 0 + 3 1 ) = 28 (2^0+2^1+2^2)\times(3^0+3^1)=28 (20+21+22)×(30+31)=28

    p 0 + p 1 + p 2 + . . . + p n = . . . ( ( 1 × p + 1 ) × p + 1 ) . . . p^0+p^1+p^2+...+p^n=...((1 \times p + 1)\times p + 1)... p0+p1+p2+...+pn=...((1×p+1)×p+1)...,其中里边迭代n次。因此可直接用结果进行迭代计算,而不用一个变量存储中间值

    最大公约数

    模板:

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

    说明:

    欧几里得算法 $\text{gcd}(a, b) = \text{gcd} (b, a\mod b) $ gcd ( a , 0 ) = a \text{gcd} (a, 0) = a gcd(a,0)=a

    欧拉函数

    φ ( n ) \varphi (n) φ(n)表示 1 1 1~ n n n中与 n n n互质的数的个数。 φ ( n ) = n ( 1 − 1 p 1 ) ( 1 − 1 p 2 ) ⋯ ( 1 − 1 p k ) \varphi \left( n \right) =n\left( 1-\frac{1}{p_1} \right) \left( 1-\frac{1}{p_2} \right) \cdots \left( 1-\frac{1}{p_k} \right) φ(n)=n(1p11)(1p21)(1pk1) 其中 n = p 1 c 1 p 2 c 2 ⋯ p k c k n=p_{1}^{c_1}p_{2}^{c_2}\cdots p_{k}^{c_k} n=p1c1p2c2pkck

    求单个欧拉函数

    模板:

    int phi(int x) { int res = x; for (int i = 2; i <= x / i; i ++ ) if (x % i == 0) { res = res / i * (i - 1); while (x % i == 0) x /= i; } if (x > 1) res = res / x * (x - 1); return res; }

    说明:

    为避免浮点数除法,编程时采用另一个形式编写代码$\varphi \left( n \right) =n\left( \frac{p_1-1}{p_1} \right) \left( \frac{p_2-1}{p_2} \right) \cdots \left( \frac{p_k-1}{p_k} \right) $可用先除后乘的方式res = res / n * (n-1)避免结果溢出例子: φ ( 6 ) = 2 \varphi(6)=2 φ(6)=2,因为 1 1 1~ 6 6 6中只有 1 1 1 5 5 5 6 6 6互质可用容斥原理证明: φ ( n ) = n − ∑ i n p i + ∑ i , j n p i p j − ∑ i , j , k n p i p j p k + ⋯ \varphi \left( n \right) =n-\sum_i{\frac{n}{p_i}+\sum_{i,j}{\frac{n}{p_ip_j}}-\sum_{i,j,k}{\frac{n}{p_ip_jp_k}}}+\cdots φ(n)=nipin+i,jpipjni,j,kpipjpkn+时间复杂度 O ( n ) O(\sqrt{n}) O(n )

    筛法求多个欧拉函数

    模板:

    int primes[N], cnt; // primes[]存储所有素数 int euler[N]; // 存储每个数的欧拉函数 bool st[N]; // st[x]存储x是否被筛掉(合数标记) void get_eulers(int n) { euler[1] = 1; for (int i = 2; i <= n; i ++ ) { if (!st[i]) { // 质数 primes[cnt ++ ] = i; euler[i] = i - 1; } for (int j = 0; primes[j] <= n / i; j ++ ) { int t = primes[j] * i; st[t] = true; if (i % primes[j] == 0) { euler[t] = euler[i] * primes[j]; break; } euler[t] = euler[i] * (primes[j] - 1); } } }

    说明:

    核心思想: φ ( n ) \varphi(n) φ(n) n n n的质指数大小无关 当 n n n是质数时, φ ( n ) = n − 1 \varphi(n)=n-1 φ(n)=n1 n m o d    p = 0 n \mod p = 0 nmodp=0时, φ ( p n ) = p φ ( n ) \varphi(pn)=p\varphi(n) φ(pn)=pφ(n),因为质数 p p p只影响 p n pn pn的质指数 c p c_p cp n m o d    p ≠ 0 n \mod p \ne 0 nmodp=0时, φ ( p n ) = p φ ( n ) ( 1 − 1 p ) = ( p − 1 ) φ ( n ) \varphi(pn)=p\varphi(n)(1-\frac{1}{p})=(p-1)\varphi(n) φ(pn)=pφ(n)(1p1)=(p1)φ(n),因为质数 p p p影响了 p n pn pn展开项数,多了一项 p 1 p^1 p1 借助线性筛法求欧拉函数,时间复杂度 O ( n ) O(n) O(n)

    快速幂

    快速计算 a k m o d    p a^k \mod p akmodp

    模板:

    typedef long long LL; int qmi(int a, int k, int p) { int res = 1; while (k) { if (k & 1) res = (LL) res * a % p; // 防止乘法溢出 a = (LL) a * a % p; // 防止乘法溢出 k >>= 1; } return res; }

    说明:

    实际上是把 k k k表示成二进制 b log k ⋯ b 2 b 1 b 0 b_{\text{log}k}\cdots b_2b_1b_0 blogkb2b1b0,因此 k = b 0 2 2 0 + b 1 2 2 1 + b 2 2 2 2 + ⋯ b log k 2 2 log k k=b_02^{2^0} + b_12^{2^1}+b_22^{2^2}+ \cdots b_{\text{log}k} 2^{2^{\text{log}k}} k=b0220+b1221+b2222+blogk22logk计算 a k m o d    p = a b 0 2 2 0 + b 1 2 2 1 + b 2 2 2 2 + ⋯ b log k 2 2 log k m o d    p = ( a b 0 2 2 0 m o d    p ) × ( a b 1 2 2 1 m o d    p ) × ⋯ × ( a b log k 2 2 log k m o d    p ) a^k \mod p=a^{b_02^{2^0} + b_12^{2^1}+b_22^{2^2}+ \cdots b_{\text{log}k} 2^{2^{\text{log}k}}} \mod p=(a^{b_0 2^{2^0}}\mod p) \times (a^{b_1 2^{2^1}}\mod p) \times \cdots \times (a^{b_{\text{log}k} 2^{2^{\text{log}k}}}\mod p) akmodp=ab0220+b1221+b2222+blogk22logkmodp=(ab0220modp)×(ab1221modp)××(ablogk22logkmodp)此外迭代满足如下关系 2 2 i + 1 m o d    p = ( 2 2 i m o d    p ) 2 2^{2^{i+1}} \mod p=(2^{2^{i}} \mod p)^2 22i+1modp=(22imodp)2中间结果可能会溢出,要强制进行类型转换时间复杂度 O ( log k ) O(\text{log}k) O(logk)

    扩展欧几里得算法

    模板

    // 求x, y,使得ax + by = gcd(a, b) int exgcd(int a, int b, int &x, int &y) { if (!b) { x = 1; y = 0; return a; } int d = exgcd(b, a % b, y, x); // 递归结束时,y = x',x = y' y -= (a/b) * x; // y = x'-(a / b) * y' = y - (a / b) * x return d; } // 详细版 // 求x, y,使得ax + by = gcd(a, b) int exgcd(int a, int b, int &x, int &y) { if (!b) { x = 1; y = 0; return a; } int x_new, y_new; int res = exgcd(b, a % b, x_new, y_new); x = y_new; // x = y' y = x_new - (a/b) * x; // y = x' - (a / b) * y' return res; }

    说明:

    在计算最大公约数时,顺便求出一组整数解 ( x , y ) (x,y) (x,y),使其满足 a x + b y = gcd ( a , b ) ax+by=\text{gcd}(a, b) ax+by=gcd(a,b)拓展:方程 a x + b y = d ax+by=d ax+by=d有解的充要条件是 gcd ( a , b ) ∣ d \text{gcd}(a, b) | d gcd(a,b)d可用于求同余方程 a x ≡ d ( m o d b ) ⇔ ∃ y ∈ Z , s . t . a x + b y = d ax\equiv d\left( \mathrm{mod}b \right) \Leftrightarrow \exists y\in \mathrm{Z},s.t. ax+by=d axd(modb)yZ,s.t.ax+by=d,但注意,用ex_gcd(a, b, x, y)求得的x不是最终x,因为此时ex_gcd求出的x是满足的 a x + b y = gcd ( a , b ) ax+by=\text{gcd}(a,b) ax+by=gcd(a,b)的,而不是满足 a x + b y = d ax+by=d ax+by=d的,二者相差 d / gcd ( a , b ) d / \text{gcd}(a, b) d/gcd(a,b)

    证明 a x + b y = g c d ( a , b ) = g c d ( b , a m o d b ) = b x ′ + a m o d b    y ′ = b x ′ + ( a − ⌊ a b ⌋ b ) y ′ = b x ′ + a y ′ − ⌊ a b ⌋ b y ′ = a y ′ + b ( x ′ − ⌊ a b ⌋ b ) \begin{aligned} ax+by&=\mathrm{gcd}\left( a,b \right) \\ &=gcd\left( b,a\mathrm{mod}b \right) \\ &=bx\prime+a\mathrm{mod}b\,\,y\prime \\ &=bx\prime+\left( a-\lfloor \frac{a}{b} \rfloor b \right) y\prime \\ &=bx\prime+ay\prime-\lfloor \frac{a}{b} \rfloor by\prime \\ &=ay\prime+b\left( x\prime-\lfloor \frac{a}{b} \rfloor b \right) \end{aligned} ax+by=gcd(a,b)=gcd(b,amodb)=bx+amodby=bx+(abab)y=bx+aybaby=ay+b(xbab) 对比 a a a b b b的系数得: { x = y ′ y = x ′ − ⌊ a b ⌋ y ′ \begin{cases} x=y\prime\\ y=x\prime-\lfloor \frac{a}{b} \rfloor y\prime\\ \end{cases} {x=yy=xbay 因此可根据递归得到的 x ′ x\prime x y ′ y\prime y计算 x x x y y y

    中国剩余定理

    已知 m 1 , m 2 ⋯ m k m_1,m_2\cdots m_k m1,m2mk互质,方程 { x ≡ a 1 ( m o d m 1 ) x ≡ a 2 ( m o d m 2 ) ⋯ x ≡ a k ( m o d m k ) \left\{ \begin{array}{c} x\equiv a_1\left( \mathrm{mod} m_1 \right)\\ \begin{array}{l} x\equiv a_2\left( \mathrm{mod} m_2 \right)\\ \cdots\\ x\equiv a_k\left( \mathrm{mod} m_k \right)\\ \end{array}\\ \end{array} \right. xa1(modm1)xa2(modm2)xak(modmk) 的最小正整数解 x = a 1 M M 1 − 1 + a 2 M M 2 − 2 + ⋯ a k M M k − k x=a_1MM_1^{-1}+a_2MM_2^{-2}+\cdots a_kMM_k^{-k} x=a1MM11+a2MM22+akMMkk,其中 M = m 1 m 2 ⋯ m k M=m_1m_2\cdots m_k M=m1m2mk M i = M m i M_i=\frac{M}{m_i} Mi=miM M i − 1 M_i^{-1} Mi1表示 M i M_i Mi m i m_i mi的逆元,可通过 M i x ≡ 1 ( m o d    m i ) M_ix\equiv1(\mod m_i) Mix1(modmi)求出

    模板

    typedef long long LL; // 扩展欧几里得算法 int exgcd(int a, int b, int &x, int &y) { if (!b) { x = 1; y = 0; return a; } int d = exgcd(b, a % b, y, x); // 递归结束时,y = x',x = y' y -= (a/b) * x; // y = x'-(a / b) * y' = y - (a / b) * x return d; } // 中国剩余定理 LL x = 0, m1, a1; cin >> m1 >> a1; for (int i = 0; i < n - 1; i ++ ) { LL m2, a2; cin >> m2 >> a2; LL k1, k2; LL d = exgcd(m1, -m2, k1, k2); // 求的不是最终解k1和k2 if ((a2 - a1) % d) { x = -1; break; // 无解 } k1 *= (a2 - a1) / d; // 变换成真正的解 int t = m2 / d; // k1 = k1 + k * m2 / d k1 = (k1 % t + t ) % t; // 变换成最小的正整数(防止C++对负数模取余) x = k1 * m1 + a1; // 把两个方程合并成一个方程 LL m = abs(m1 / d * m2); // 变成正数 a1 = k1 * m1 + a1; m1 = m; } if (x != -1) x = (x % m1 + m1) % m1; // 变成最小正整数(防止C++对负数模取余) cout << x << endl;

    说明

    本质是求解同余方程组如果直接按公式求解,结果可能会溢出,因此采用迭代的方式求解 x ≡ a ( m o d m ) ⇔ x = m k + a , k ∈ Z x\equiv a\left( \mathrm{mod} m \right) \Leftrightarrow x=mk+a,k\in \mathrm{Z} xa(modm)x=mk+a,kZ考虑第1个和第2个方程,可得等式 m 1 k 1 + a 1 = m 2 k 2 + a 2 m_1k_1+a_1=m_2k_2+a_2 m1k1+a1=m2k2+a2,变形得 m 1 k 1 − m 2 k 2 = a 2 − a 1 m_1k_1-m_2k_2=a_2-a_1 m1k1m2k2=a2a1 d = gcd ( m 1 , − m 2 ) d=\text{gcd}(m_1, -m_2) d=gcd(m1,m2)方程有解的充分必要条件是 d ∣ ( a 2 − a 1 ) d | (a_2-a_1) d(a2a1)不定方程的一个通解是 { k 1 = k 1 + k m 2 d k 2 = k 2 + k m 1 d \begin{cases} k_1=k_1+k\frac{m_2}{d}\\ k_2=k_2+k\frac{m_1}{d}\\ \end{cases} {k1=k1+kdm2k2=k2+kdm1,其中 k k k是参数此时 x = k 1 m 1 + a 1 = ( k 1 + k m 2 d ) m 1 + a 1 = k 1 m 1 + k m 1 m 2 d + a 1 = k 1 m 1 + a 1 + k m 1 m 2 d = x 1 + k d 12 x=k_1m_1+a_1=(k_1+k\frac{m_2}{d})m_1+a_1=k_1m_1+k\frac{m_1m_2}{d}+a_1=k_1m_1+a_1+k\frac{m_1m_2}{d}=x_1+kd_{12} x=k1m1+a1=(k1+kdm2)m1+a1=k1m1+kdm1m2+a1=k1m1+a1+kdm1m2=x1+kd12,这里令 x 1 = k 1 m 1 + a 1 x_1=k_1m_1+a_1 x1=k1m1+a1 d 12 = m 1 m 2 d d_{12}=\frac{m_1m_2}{d} d12=dm1m2若把 x 1 x_1 x1看成 a 1 a_1 a1 d 12 d_{12} d12看成 m 1 m_1 m1,则上式变成 x = a 1 + k m 1 x=a_1+km_1 x=a1+km1,即 x ≡ a 1 ( m o d    m 1 ) x\equiv a_1 (\mod m_1) xa1(modm1),相当于把两个同余方程合并成了一个

    高斯消元

    模板

    // a[N][N+1]是增广矩阵 int gauss() { int c, r; // 按列遍历,化成行阶梯矩阵 for (c = 0, r = 0; c < n; c ++) { // 找到当前列绝对值最大元素所在的行(搜索第r行~最后一行) int t = r; for (int i = r + 1; i < n; i ++ ) if (fabs(a[i][c]) > fabs(a[t][c])) t = i; if (fabs(a[t][c]) < eps) continue; // 当前列全是0 for (int j = c; j <= n; j ++ ) swap(a[t][j], a[r][j]); // 将绝对值最大的行换到最顶端 for (int j = n; j >= c; j -- ) a[r][j] /= a[r][c]; // 将当前上的首位变成1,注意从后往前遍历 for (int i = r + 1; i < n; i ++ ) // 用当前行将下面所有的列消成0 if (fabs(a[i][c]) > eps) for (int j = n; j >= c; j -- ) a[i][j] -= a[r][j] * a[i][c]; // 注意从后往前删,否则出现写后读错误 r ++ ; } if (r < n) { for (int i = r; i < n; i ++ ) if (fabs(a[i][n]) > eps) return 2; // 出现0!=0,无解 return 1; // 都是0=0,有无穷多组解 } // 化成单位阵,增广矩阵的扩展部分为方程的解 for (int i = n - 1; i >= 0; i -- ) for (int j = i + 1; j < n; j ++ ) a[i][n] -= a[i][j] * a[j][n]; return 0; // 有唯一解 }

    说明

    可在 O ( n 3 ) O(n^3) O(n3)求解线性方程组 找到绝对值最大的一行将该行换到最上面将该行第一个数化成 1 1 1把该列下边的数化成 0 0 0 注意部分过程需要从后往前遍历,否则出现写后读问题

    组合数

    递归法

    C a b = C a − 1 b + C a − 1 b − 1 C_{a}^{b}=C_{a-1}^{b}+C_{a-1}^{b-1} Cab=Ca1b+Ca1b1

    模板

    for (int i = 0; i < n; i ++ ) for (int j = 0; j <= i; j ++ ) if (!j) c[i][j] = 1; else c[i][j] = c[i - 1][j] + c[i - 1][j - 1];

    说明

    意义: C a b C_{a}^{b} Cab表示从 a a a个苹果中选 b b b个的方案数i=0时,不会执行else部分,因此不会出现数组越界代码结构有点像杨辉三角适用于询问次数多,但组合数取值范围小的情况( 1 0 3 10^3 103数量级)

    预处理逆元

    假设 m m m是一个很大的质数

    ( n − 1 ) ! (n-1)! (n1)!的模 m m m逆元是 ( ( n − 1 ) ! ) m − 2 ((n-1)!)^{m-2} ((n1)!)m2

    n ! n! n!的模 m m m逆元是 ( n ! ) m − 2 (n!)^{m-2} (n!)m2

    infact ( n ) = infact ( n − 1 ) × ( n ! ) m − 2 ( ( n − 1 ) ! ) m − 2 = infact ( n − 1 ) × n m − 2 \text{infact}(n)=\text{infact}(n-1) \times \frac{(n!)^{m-2}}{((n-1)!)^{m-2}}=\text{infact}(n-1) \times n^{m-2} infact(n)=infact(n1)×((n1)!)m2(n!)m2=infact(n1)×nm2

    其中 n m − 2 n^{m-2} nm2可通过快速幂求得

    模板

    // 快速幂模板 int qmi(int a, int k, int p) { int res = 1; while (k) { if (k & 1) res = (LL)res * a % p; a = (LL)a * a % p; k >>= 1; } return res; } // 预处理阶乘的余数和阶乘逆元的余数 fact[0] = infact[0] = 1; for (int i = 1; i < N; i ++ ) { fact[i] = (LL) fact[i - 1] * i % mod; infact[i] = (LL) infact[i - 1] * qmi(i, mod - 2, mod) % mod; }

    说明

    在涉及乘法运算时,注意先强制转换成long long,然后再取模,否则结果会溢出适用于询问次数多,但组合数取值范围较大的情况( 1 0 6 10^6 106数量级)

    Lucas定理

    p p p是质数,则对于任意整数 1 ≤ m ≤ n 1 \le m \le n 1mn,有: C n m = C n m o d    p m m o d    p × C n ÷ p m ÷ p m o d    p C_n^m=C_{n \mod p}^{m \mod p} \times C_{n \div p}^{m \div p} \mod p Cnm=Cnmodpmmodp×Cn÷pm÷pmodp

    模板

    // 快速幂模板 int qmi(int a, int k) { int res = 1; while (k) { if (k & 1) res = (LL)res * a % p; a = (LL)a * a % p; k >>= 1; } return res; } // 通过定理求组合数C(a, b) int C(int a, int b) { int res = 1; for (int i = 1, j = a; i <= b; i ++, j -- ) { res = (LL)res * j % p; // 构造分子a * (a - 1) * ... * (a - b + 1) res = (LL)res * qmi(i, p - 2) % p; // 构造(b!)的逆元(b!)^{p-2} } return res; } int lucas(LL a, LL b) { if (a < p && b < p) return C(a, b); return (LL)C(a % p, b % p) * lucas(a / p, b / p) % p; }

    说明

    适用于询问次数少,但组合数取值范围较大的情况( a , b ≤ 1 0 18 a, b \le10^{18} a,b1018 p ≤ 1 0 5 p\le 10^5 p105)两方面优化组合数的求解 快速幂加速计算小规模组合数卢卡斯定理降低组合数的规模 注意变量类型以及输入输出的类型是int还是long long

    分解质因数法

    模板

    int primes[N], cnt; // 存储所有质数 int sum[N]; // 存储每个质数的次数 bool st[N]; // 存储每个数是否已被筛掉(合数标记) // 线性筛法求素数 void get_primes(int n) { for (int i = 2; i <= n; i ++ ) { if (!st[i]) primes[cnt ++ ] = i; for (int j = 0; primes[j] <= n / i; j ++ ) { st[primes[j] * i] = true; if (i % primes[j] == 0) break; } } } // 求n!质因数分解后,在质数p的次数 int get(int n, int p) { int res = 0; while (n) { res += n / p; n /= p; } return res; } // 高精度乘低精度模板 vector<int> mul(vector<int> a, int b) { vector<int> c; int t = 0; for (int i = 0; i < a.size(); i ++ ) { t += a[i] * b; c.push_back(t % 10); t /= 10; } while (t) { c.push_back(t % 10); t /= 10; } return c; } get_primes(a); // 预处理范围内的所有质数 for (int i = 0; i < cnt; i ++ ) // 求每个质因数的次数 { int p = primes[i]; sum[i] = get(a, p) - get(b, p) - get(a - b, p); } vector<int> res; res.push_back(1); for (int i = 0; i < cnt; i ++ ) // 用高精度乘法将所有质因子相乘 for (int j = 0; j < sum[i]; j ++ ) res = mul(res, primes[i]);

    说明

    本质是质因数分解+高精度乘法适用求真实值,而非取余后的结果主要思想 线性筛法找出 2 2 2~ a a a内的所有质数依次遍历质数,求出组合数在各个质数的次数(通过get方法计算)高精度乘法计算各个质因子的乘积 get(int n, int p)求n!质因数分解在 p p p的次数的数学公式$\lfloor \frac{n}{p} \rfloor +\lfloor \frac{n}{p^2} \rfloor +\lfloor \frac{n}{p^3} \rfloor +\cdots $

    卡特兰数

    catalan ( n ) = C 2 n n n + 1 \text{catalan}(n)=\frac{C_{2n}^n}{n+1} catalan(n)=n+1C2nn

    容斥原理

    ∣ ⋃ k = 1 n S k ∣ = ∑ i ∣ S i ∣ − ∑ i , j ∣ S i ∩ S j ∣ + ∑ i , j , k ∣ S i ∩ S j ∩ S k ∣ − ⋯ \left| \bigcup_{k=1}^n{S_k} \right|=\sum_i{\left| S_i \right|}-\sum_{i,j}{\left| S_i\cap S_j \right|}+\sum_{i,j,k}{\left| S_i\cap S_j\cap S_k \right|}-\cdots k=1nSk=iSii,jSiSj+i,j,kSiSjSk

    C n 1 + C n 2 + C n 3 + ⋯ + C n n = 2 n − C n 0 = 2 n − 1 C_{n}^{1}+C_{n}^{2}+C_{n}^{3}+\cdots +C_{n}^{n}=2^n-C_{n}^{0}=2^n-1 Cn1+Cn2+Cn3++Cnn=2nCn0=2n1

    模板

    int res = 0; // 一共有2^m-1种集合 for (int i = 1; i < 1 << m; i ++ ) { int t = 1, s = 0; // t表示当前集合的乘积,s表示符号(正负交替) // 遍历每种集合 for (int j = 0; j < m; j ++ ) if (i >> j & 1) { // 元素存在 if ((LL)t * p[j] > n) { t = -1; // 构造的质数乘积比n大,舍弃 break; } t *= p[j]; s ++ ; } if (t != -1) { if (s % 2) res += n / t; // 根据s的奇偶性实现正负交替 else res -= n / t; } }

    说明

    2 m 2^m 2m可以用1 << m实现

    m个元素可以构造 2 m − 1 2^m-1 2m1种不同的非空集合,集合包含的元素可用 m m m位二进制表示,例如5=101b,表示集合有第1个元素和第3个元素

    对于每个元素都有 C n 1 − C n 2 + C n 3 − ⋯ + ( − 1 ) k − 1 C n n = 1 C_{n}^{1}-C_{n}^{2}+C_{n}^{3}-\cdots +\left( -1 \right) ^{k-1}C_{n}^{n}=1 Cn1Cn2+Cn3+(1)k1Cnn=1,即每个元素取到的次数都是1,不重复

    1 1 1~ n n n中能被 p p p整除的数的个数为$\lfloor \frac{n}{p} \rfloor $

    遍历的次序与公式不太一样,不是先加 m m m个数,再减去 C m 2 C_m^2 Cm2个数。实际上,当集合元素个数是奇数时,符号为正,反之为负。因此可根据集合个数的奇偶性判断符号的正负

    模板是基于整除个数问题设计的

    时间复杂度为 O ( 2 n ) O(2^n) O(2n)

    NIM游戏

    公平组合游戏ICG 若一个游戏满足以下条件,则称该游戏为一个公平组合游戏。

    由两名玩家交替行动;在游戏进程的任意时刻,可以执行的合法行动与轮到哪名玩家无关;不能行动的玩家判负;

    NIM博弈属于公平组合游戏,但城建的棋类游戏,比如围棋,就不是公平组合游戏。因为围棋交战双方分别只能落黑子和白子,胜负判定也比较复杂,不满足条件2和条件3。

    取石子

    给定N堆物品,第 i i i堆物品有 A i A_i Ai个。两名玩家轮流行动,每次可以任选一堆,取走任意多个物品,可把一堆取光,但不能不取。取走最后一件物品者获胜。两人都采取最优策略,问先手是否必胜。

    定义

    我们把这种游戏称为NIM博弈。把游戏过程中面临的状态称为局面。整局游戏第一个行动的称为先手,第二个行动的称为后手。若在某一局面下无论采取何种行动,都会输掉游戏,则称该局面必败。 所谓采取最优策略是指,若在某一局面下存在某种行动,使得行动后对面面临必败局面,则优先采取该行动。同时,这样的局面被称为必胜。我们讨论的博弈问题一般都只考虑理想情况,即两人均无失误,都采取最优策略行动时游戏的结果。 NIM博弈不存在平局,只有先手必胜和先手必败两种情况。

    必胜状态,先手进行某一个操作,留给后手是一个必败状态时,对于先手来说是一个必胜状态。即先手可以走到某一个必败状态。 必败状态,先手无论如何操作,留给后手都是一个必胜状态时,对于先手来说是一个必败状态。即先手走不到任何一个必败状态。

    定理

    NIM博弈先手必胜,当且仅当 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A n ≠ 0 A_1\oplus A_2\oplus \cdots \oplus A_n\ne 0 A1A2An=0

    性质

    操作到最后时,每堆石子数都是0, 0 ⊕ 0 ⊕ ⋯ ⊕ 0 = 0 0\oplus 0\oplus \cdots \oplus 0=0 000=0在操作过程中,如果 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A n = x ≠ 0 A_1\oplus A_2\oplus \cdots \oplus A_n=x\ne 0 A1A2An=x=0。那么玩家必然可以通过拿走某一堆若干个石子将异或结果变为0。在操作过程中,如果 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A n = 0 A_1\oplus A_2\oplus \cdots \oplus A_n= 0 A1A2An=0,那么无论玩家怎么拿,必然会导致最终异或结果不为 0 0 0

    性质证明

    性质2的证明:不妨设 x x x的二进制表示中最高一位 1 1 1在第 k k k位,那么在 A 1 , A 2 , ⋯   , A n A_1,A_2,\cdots,A_n A1,A2,,An中,必然有一个数 A i A_i Ai,它的第 k k k为时 1 1 1,且 A i ⊕ x < A i A_i\oplus x<A_i Aix<Ai,那么从第 i i i堆石子中拿走 A i − A i ⊕ x A_i-A_i\oplus x AiAix个石子,第ii堆石子还剩$A_i-\left( A_i-A_i\oplus x \right) , 此 时 ,此时 A_1\oplus A_2\oplus \cdots \oplus A_i\oplus x\oplus \cdots \oplus A_n=x\oplus x=0$性质3的证明:反证法:假设玩家从第 i i i堆石子拿走若干个,结果仍是 0 0 0。不妨设还剩下 A i ′ A_i′ Ai个,因为不能不拿,所以 0 ≤ A i ′ < A i 0≤A_i′<A_i 0Ai<Ai,且 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A i ′ ⊕ ⋯ ⊕ A n = 0 A_1\oplus A_2\oplus \cdots \oplus A_i′\oplus \cdots \oplus A_n=0 A1A2AiAn=0。那么 ( A 1 ⊕ A 2 ⊕ ⋯ ⊕ A i ⊕ ⋯ ⊕ A n ) ⊕ ( A 1 ⊕ A 2 ⊕ ⋯ ⊕ A i ′ ⊕ ⋯ ⊕ A n ) = A i ⊕ A i ′ = 0 \left( A_1\oplus A_2\oplus \cdots \oplus A_i\oplus \cdots \oplus A_n \right) \oplus \left( A_1\oplus A_2\oplus \cdots \oplus A_i′\oplus \cdots \oplus A_n \right) =A_i\oplus A_i′=0 (A1A2AiAn)(A1A2AiAn)=AiAi=0,则 A i = A i ′ A_i= A_i′ Ai=Ai,与假设 0 ≤ A i ′ < A i 0≤A_i′<A_i 0Ai<Ai矛盾。

    分析

    当先手面对的局面是 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A n ≠ 0 A_1\oplus A_2\oplus \cdots \oplus A_n\ne 0 A1A2An=0,则先手总有办法让后手局面变成 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A n = 0 A_1\oplus A_2\oplus \cdots \oplus A_n= 0 A1A2An=0,后手必败,先手必赢当先手面对的局面是 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A n = 0 A_1\oplus A_2\oplus \cdots \oplus A_n= 0 A1A2An=0,则后手总有办法让先手局面变成 A 1 ⊕ A 2 ⊕ ⋯ ⊕ A n ≠ 0 A_1\oplus A_2\oplus \cdots \oplus A_n\ne 0 A1A2An=0,先手必败,后手必赢

    取石子(升级版)

    约束条件

    每次可拿石子的个数不是任意的,而是集合 S S S中的元素

    概念

    mex \text{mex} mex

    mex(S) \text{mex(S)} mex(S)表示不属于集合 S S S的最小自然数,即$\text{mex}(S)=\min \left{ x \mid x\notin S\land x\in \mathrm{N} \right} $。

    例如当 S = 1 , 2 S={1, 2} S=1,2时, mex ( S ) = 0 \text{mex}(S)=0 mex(S)=0例如当 S = 0 , 2 S={0, 2} S=0,2时, mex ( S ) = 1 \text{mex}(S)=1 mex(S)=1

    有向图游戏

    给定一个有向无环图,图中有一个唯一的起点,在起点上放有一枚棋子。两名玩家交替地把这枚棋子沿有向边进行移动,每次可以移动一步,无法移动者判负。该游戏被称为有向图游戏。

    任何一个公平组合游戏都可以转化为有向图游戏。具体方法是,把每个局面看成图中的一个节点,并且从每个局面向沿着合法行动能够到达的下一个局面连有向边。

    SG \text{SG} SG

    在有向图游戏中,对于每个节点 x x x,设从 x x x出发共有 k k k条有向边,分别到达节点 y 1 , y 2 , ⋯   , y k y_1, y_2,\cdots, y_k y1,y2,,yk,定义 SG ( x ) \text{SG}(x) SG(x) x x x的后继节点 y 1 , y 2 , ⋯   , y k y_1, y_2,\cdots, y_k y1,y2,,yk SG \text{SG} SG函数值构成的集合再执行 mex ( S ) \text{mex}(S) mex(S)运算的结果,即: SG ( x ) = mex ( SG ( y 1 ) , SG ( y 2 ) , ⋯   , SG ( y k ) ) \text{SG}(x) = \text{mex}({\text{SG}(y1), \text{SG}(y2),\cdots, \text{SG}(yk)}) SG(x)=mex(SG(y1),SG(y2),,SG(yk)) 特别地,整个有向图游戏 G G G SG \text{SG} SG函数值被定义为有向图游戏起点 s s s SG \text{SG} SG函数值,即 SG ( G ) = SG ( s ) \text{SG}(G) = \text{SG}(s) SG(G)=SG(s)

    有向图游戏的和

    G 1 , G 2 , ⋯   , G m G_1, G_2, \cdots, G_m G1,G2,,Gm m m m个有向图游戏。定义有向图游戏 G G G,它的行动规则是任选某个有向图游戏 G i G_i Gi,并在 G i G_i Gi上行动一步。 G G G被称为有向图游戏 G 1 , G 2 , ⋯   , G m G_1, G_2,\cdots, G_m G1,G2,,Gm的和。 有向图游戏的和的SG函数值等于它包含的各个子游戏SG函数值的异或和,即: S G ( G ) = S G ( G 1 ) ⊕ S G ( G 2 ) ⊕ ⋯ ⊕ S G ( G m ) \mathrm{SG}\left( G \right) =\mathrm{SG}\left( G_1 \right) \oplus \mathrm{SG}\left( G_2 \right) \oplus \cdots \oplus \mathrm{SG}\left( G_m \right) SG(G)=SG(G1)SG(G2)SG(Gm)

    性质

    SG ( k ) 有 \text{SG}(k)有 SG(k) k k k个后继结点,分别是 0 0 0~ k − 1 k-1 k1 0 0 0可以走向 0 0 0 0 0 0只能走向非 0 0 0

    性质证明

    类似NIM游戏的证明

    定理

    有向图游戏的某个局面必胜,当且仅当该局面对应节点的SG函数值大于 0 0 0。 有向图游戏的某个局面必败,当且仅当该局面对应节点的SG函数值等于 0 0 0

    模板

    int s[N], f[M]; memset(f, -1, sizeof f); // 初始化 // 记忆化搜索(备忘录法) int sg(int x) { if (f[x] != -1) return f[x]; // 已经计算过 // 构造子树 unordered_set<int> S; for (int i = 0; i < m; i ++ ) { int sum = s[i]; if (x >= sum) S.insert(sg(x - sum)); // 保存后继结点到S中(递归) } // mem(x) for (int i = 0; ; i ++ ) if (!S.count(i)) return f[x] = i; } // 把n堆石子看成n个独立的有向图G,把各个有向图结果做异或即可得到答案 int x, res = 0; for (int i = 0; i < n; i ++ ) { cin >> x; res ^= sg(x); }

    P.S. 部分内容来自y总的模板 如果大家有兴趣,可以去Acwing《算法基础课》看看 我在Acwing也分享了一份,欢迎去围观

    Processed: 0.015, SQL: 8