数论——同余问题全家桶3 __int128和同余方程组
数论——同余问题全家桶3 __int128和同余方程组
- 快速读写和__int128
- 快速读写
- `__int128`
- 中国剩余定理和线性同余方程组
- 中国剩余定理(CRT)
- 中国剩余定理OJ示例
- 模板题曹冲养猪 - 洛谷
- 模板题猜数字 - 洛谷
- 扩展中国剩余定理
- 扩展中国剩余定理OJ示例
- 模板题扩展中国剩余定理(EXCRT) - 洛谷
- NOI2018 屠龙勇士
快速读写和__int128
这块算是补充知识点,但仅作为高精度算法的临时替代,在有的题若使用
__int128
也会溢出,则只能用高精度算法。
虽然在 99%
的题目中,scanf
以及 printf
来处理输入和输出已经足够快了。但是,如果输入和输出的规模特别庞大,或者出题人卡常,scanf
和 printf
也是会超时的,这个时候就需要更加快速的读写方式。
同样的,虽然在 99%
的题目中,long long
已经足够我们应付较大的整数。但是,如果题目中数的范围过大,相乘也是有可能超过 long long
的最大范围。如果只是大一点点,此时可以用 __int128
来存储。
快速读写
快速读写有很多种版本,接下来要介绍一种容易实现且够用的版本 - getchar
/putchar
结合秦九韶算法。
前置知识:
-
在计算机的视角,所有的整数其实是一个一个的字符串,每个数拆开看就是一个一个字符。因此,对于一个整数,可以当成字符串,一个一个字符的输入进来。同理,输出一个整数的时候,也可以当成字符串,一个一个字符的输出。
-
getchar
/putchar
这种输入输出方式相较于scanf
/printf
而言,速度更快。
于是,就可以利用 getchar
将字符转换成整数输入,利用 putchar
将整数转换成字符输出。
如果这种方式还是超时,那就把 getchar
换成 getchar_unlocked
。如果换完之后还是超时,那就是毒瘤题,可以忽视,或者就是算法本身的时间复杂度有问题。
但是,getchar_unlocked
只能在 Linux 系统下使用,因此要慎用。
这里通过模板题进行检测:
P10815 【模板】快速读入 - 洛谷
#include<bits/stdc++.h>
using namespace std;inline int read(){int flag=1,ans=0;//这个函数只有在linux操作系统下的gcc系列编译器才能使用 //当然,洛谷也能使用,不然最后一个样例无法通过 char ch=getchar_unlocked();//一般编译器用这个
// char ch=getchar();while(ch<'0'||ch>'9'){if(ch=='-')flag=-1;ch=getchar_unlocked();
// ch=getchar();}while(ch>='0'&&ch<='9'){ans=ans*10+ch-'0';ch=getchar_unlocked();
// ch=getchar();}return ans*flag;
}inline void print(int a){if(a<0){putchar('-');a=-a;}if(a>9)print(a/10);putchar(a%10+'0');
}int main() {int n=read();int sum=0,x;while(n--){x=read();sum+=x;}print(sum);return 0;
}
__int128
【更大的整数__int128】
__int128
就是占用 128 字节的整数存储类型,范围就是 − 2 127 ∼ 2 127 − 1 -2^{127} \sim 2^{127} - 1 −2127∼2127−1。
即
− 170141183460469231731687303715884105728 ∼ 170141183460469231731687303715884105727 -170141183460469231731687303715884105728\\\sim170141183460469231731687303715884105727 −170141183460469231731687303715884105728∼170141183460469231731687303715884105727。
170 1411 8346 0469 2317 3168 7303 7158 8410 5727//39位
-
对比
int
的范围 − 2 32 ∼ 2 32 − 1 -2^{32} \sim 2^{32} - 1 −232∼232−1,即 − 2147483648 ∼ 2147483647 -2147483648\sim2147483647 −2147483648∼2147483647。无符号情况下 0 ∼ 4294967295 0\sim 4294967295 0∼4294967295。
-
long long
: − 2 64 ∼ 2 64 − 1 -2^{64} \sim 2^{64} - 1 −264∼264−1,即 − 9223372036854775808 ∼ 9223372036854775807 -9223372036854775808\sim9223372036854775807 −9223372036854775808∼9223372036854775807。无符号情况下 0 ∼ 18446744073709551615 0\sim18446744073709551615 0∼18446744073709551615。
- 如果使用了
unsigned __int128
,则范围变成 0 ∼ 2 128 0 \sim 2^{128} 0∼2128,即约 39 位数。
0 ∼ 340282366920938463463374607431768211455 0\sim 340282366920938463463374607431768211455 0∼340282366920938463463374607431768211455。
当数据范围超过 long long
,但是还不足以用上高精度算法的时候,用 __int128
是个不错的选择。
但是,__int128
是不能直接用 cin
、cout
、scanf
以及 printf
输入或者输出的。只能按照字符的方式输入输出,也就是我们刚刚学习的快速读写的方式。不过,当读入进来之后,用法和普通的整型变量一致。
__int128
的使用示例(需要在gcc系IDE使用,例如Devc++):
#include<bits/stdc++.h>
using namespace std;
typedef unsigned __int128 uint;
typedef __int128 _int;inline void print(uint n){if(n>9)print(n/10);putchar(n%10+'0');
}inline void print(_int n){if(n<0){putchar('-');n=-n;}if(n>9)print(n/10);putchar(n%10+'0');
}inline _int read(){_int flag=1,ans=0;char ch=getchar();while(ch<'0'||ch>'9'){if(ch=='-')flag=-1;ch=getchar();}while(ch>='0'&&ch<='9'){ans=ans*10+ch-'0';ch=getchar();}return ans*flag;
}void f1(){_int mmax=(_int(1)<<127)-1;print(mmax);cout<<endl;_int mmin=(_int(1)<<127);print(mmin+1);//无法直接输出-2^128 cout<<endl;
}void f2(){uint mmax=~uint(0);print(mmax);cout<<endl;
}void f3(){_int a=114514;print(a);cout<<endl;a=read();print(a);cout<<endl;a+=3;print(a);cout<<endl;a-=486;print(a);cout<<endl;a*=1435;print(a);cout<<endl;a/=1435;print(a);cout<<endl;a%=10007;print(a);cout<<endl;
}int main() {
// f1();
// f2();f3();return 0;
}
__int128
再好用,缺陷也很多:
-
不通用。
__int128
并没有在任何一个 c++ 标准中严格定义,所以目前它只是 GCC 系列编译器的专属。目前测试,只在 Linux 系统下能够正常使用,在其他编译器例如MSVC不能使用。因此如果要使用,就要看比赛中的编译器,是否是 Linux。 -
不方便。
__int128
目前是不支持直接读入、输出的。也就是无法用cin
、cout
、scanf
以及printf
输入或者输出这种类型的数。只能按照字符的方式输入输出,也就是我们刚刚学习的快速读写的方式。 -
空间大,速度慢。
__int128
占用了 16 个字节来存,空间超限的风险显著增加。
但是,在有些题目中,__int128
还是足够我们使用的。
中国剩余定理和线性同余方程组
引入线性同余方程组:南北朝时期《孙子算经》中有个问题:有物不知其数,三三数之剩二,五五数之剩三,七七数之剩二。问物几何?
用数学公式翻译过来就是一个线性同余方程组:
{ x ≡ 2 ( mod 3 ) x ≡ 3 ( mod 5 ) x ≡ 2 ( mod 7 ) \begin{cases} x \equiv 2(\text{mod } 3) \\ x \equiv 3(\text{mod } 5) \\ x \equiv 2(\text{mod } 7) \end{cases} ⎩ ⎨ ⎧x≡2(mod 3)x≡3(mod 5)x≡2(mod 7)
线性同余方程组:求这个方程组的最小非负整数解。
{ x ≡ r 1 ( mod m 1 ) x ≡ r 2 ( mod m 2 ) ⋮ x ≡ r n ( mod m n ) \begin{cases} x \equiv r_1(\text{mod } m_1) \\ x \equiv r_2(\text{mod } m_2) \\ \vdots \\ x \equiv r_n(\text{mod } m_n) \end{cases} ⎩ ⎨ ⎧x≡r1(mod m1)x≡r2(mod m2)⋮x≡rn(mod mn)
中国剩余定理(CRT)
前提:所有的模数 m 1 , m 2 , . . . , m n m_1, m_2, ..., m_n m1,m2,...,mn 两两互质。因此中国剩余定理的使用比较局限。
原理:中国剩余定理是基于“构造法”得出的结果。以下给出 n = 3 n = 3 n=3 的构造过程, n n n 等于任意数的构造过程是一样的。
对于这个方程组:
{ x ≡ r 1 ( mod m 1 ) x ≡ r 2 ( mod m 2 ) x ≡ r 3 ( mod m 3 ) \begin{cases} x \equiv r_1 (\text{mod } m_1) \\ x \equiv r_2 (\text{mod } m_2)\\x \equiv r_3 (\text{mod } m_3) \end{cases} ⎩ ⎨ ⎧x≡r1(mod m1)x≡r2(mod m2)x≡r3(mod m3)
或
{ x m o d m 1 = r 1 x m o d m 2 = r 2 x m o d m 3 = r 3 \begin{cases}x\bmod m_1=r_1\\x\bmod m_2=r_2\\x\bmod m_3=r_3\end{cases} ⎩ ⎨ ⎧xmodm1=r1xmodm2=r2xmodm3=r3
记 M = m 1 × m 2 × m 3 M = m_1 \times m_2 \times m_3 M=m1×m2×m3,构造解: x = C 1 + C 2 + C 3 x = C_1 + C_2 + C_3 x=C1+C2+C3,其中:
-
C 1 mod m 1 = r 1 , C 1 mod m 2 = 0 , C 1 mod m 3 = 0 C_1 \text{ mod } m_1 = r_1,\ C_1 \text{ mod } m_2 = 0,\ C_1 \text{ mod } m_3 = 0 C1 mod m1=r1, C1 mod m2=0, C1 mod m3=0
-
C 2 mod m 1 = 0 , C 2 mod m 2 = r 2 , C 2 mod m 3 = 0 C_2 \text{ mod } m_1 = 0, \ C_2 \text{ mod } m_2 = r_2, \ C_2 \text{ mod } m_3 = 0 C2 mod m1=0, C2 mod m2=r2, C2 mod m3=0
-
C 3 mod m 1 = 0 , C 3 mod m 2 = 0 , C 3 mod m 3 = r 3 C_3 \text{ mod } m_1 = 0, \ C_3 \text{ mod } m_2 = 0, \ C_3 \text{ mod } m_3 = r_3 C3 mod m1=0, C3 mod m2=0, C3 mod m3=r3
这样 x m o d m 1 = C 1 m o d m 1 + C 2 m o d m 1 + C 3 m o d m 1 = r 1 x\bmod m1 = C_1\bmod m1+C_2\bmod m1+C_3\bmod m_1=r_1 xmodm1=C1modm1+C2modm1+C3modm1=r1,其他 m 2 m2 m2、 m 3 m3 m3同理。
因此只要能构造出这样的 x x x,那么 x x x 就是我们要的解。
C 1 C_1 C1、 C 2 C_2 C2和 C 3 C_3 C3的构造方式:
-
C 1 = r 1 × m 2 × m 3 × ( m 2 × m 3 ) − 1 , C_1 = r_1 \times m_2 \times m_3 \times (m_2 \times m_3)^{-1}, C1=r1×m2×m3×(m2×m3)−1,
$(m_2 \times m_3)^{-1} 为模 为模 为模m_1 意义下的逆元, 意义下的逆元, 意义下的逆元,m_2 \times m_3 \times (m_2 \times m_3)^{-1}\bmod m_1 最后的结果是 1 。所以 最后的结果是1。所以 最后的结果是1。所以C_1\bmod m_1=r_1 , , ,C_1\bmod m_2=C_1\bmod m_3=0$。
-
C 2 = r 2 × m 1 × m 3 × ( m 1 × m 3 ) − 1 , C_2 = r_2 \times m_1 \times m_3 \times (m_1 \times m_3)^{-1}, C2=r2×m1×m3×(m1×m3)−1, 其中 $(m_1 \times m_3)^{-1} 为模 为模 为模m_1$意义下的逆元。
-
C 3 = r 3 × m 1 × m 2 × ( m 1 × m 2 ) − 1 , C_3 = r_3 \times m_1 \times m_2 \times (m_1 \times m_2)^{-1}, C3=r3×m1×m2×(m1×m2)−1, 其中 $(m_1 \times m_2)^{-1} 为模 为模 为模m_1$意义下的逆元。
当 x x x加减 M M M的若干倍时,依旧是满足方程,因此最小非负整数解就是
( x mod M + M ) mod M (x \text{ mod } M + M) \text{ mod } M (x mod M+M) mod M。
整个方程的通解是 x + k M x+kM x+kM。
以《孙子算经》中的问题为例:
{ x ≡ 2 ( mod 3 ) x ≡ 3 ( mod 5 ) x ≡ 2 ( mod 7 ) \begin{cases} x \equiv 2(\text{mod } 3) \\ x \equiv 3(\text{mod } 5) \\ x \equiv 2(\text{mod } 7) \end{cases} ⎩ ⎨ ⎧x≡2(mod 3)x≡3(mod 5)x≡2(mod 7)
-
计算每一个方程的 C i C_i Ci:
C 1 = r 1 × m 2 × m 3 × ( m 2 × m 3 ) − 1 = 2 × 5 × 7 × 35 − 1 ( mod 3 ) = 2 × 5 × 7 × 2 = 140 C_1 = r_1 \times m_2 \times m_3 \times (m_2 \times m_3)^{-1} = 2 \times 5 \times 7 \times 35^{-1}(\text{mod } 3) = 2 \times 5 \times 7 \times 2 = 140 C1=r1×m2×m3×(m2×m3)−1=2×5×7×35−1(mod 3)=2×5×7×2=140
C 2 = r 2 × m 1 × m 3 × ( m 1 × m 3 ) − 1 = 3 × 3 × 7 × 21 − 1 ( mod 5 ) = 3 × 3 × 7 × 1 = 63 C_2 = r_2 \times m_1 \times m_3 \times (m_1 \times m_3)^{-1} = 3 \times 3 \times 7 \times 21^{-1}(\text{mod } 5) = 3 \times 3 \times 7 \times 1 = 63 C2=r2×m1×m3×(m1×m3)−1=3×3×7×21−1(mod 5)=3×3×7×1=63
C 3 = r 3 × m 1 × m 2 × ( m 1 × m 2 ) − 1 = 2 × 3 × 5 × 15 − 1 ( mod 7 ) = 2 × 3 × 5 × 1 = 30 C_3 = r_3 \times m_1 \times m_2 \times (m_1 \times m_2)^{-1} = 2 \times 3 \times 5 \times 15^{-1}(\text{mod } 7) = 2 \times 3 \times 5 \times 1 = 30 C3=r3×m1×m2×(m1×m2)−1=2×3×5×15−1(mod 7)=2×3×5×1=30 -
计算结果: x = ( 140 + 63 + 30 ) mod 105 = 233 mod 105 = 23 x = (140 + 63 + 30) \text{ mod } 105 = 233 \text{ mod } 105 = 23 x=(140+63+30) mod 105=233 mod 105=23
推广: n n n取任意数的时候,构造方式也是同理。
- C 1 = r 1 × m 2 × m 3 × … × m n × ( m 2 × m 3 × … × m n ) − 1 C_1 = r_1 \times m_2 \times m_3 \times \ldots \times m_n \times (m_2 \times m_3 \times \ldots \times m_n)^{-1} C1=r1×m2×m3×…×mn×(m2×m3×…×mn)−1,
- C 2 = r 2 × m 1 × m 3 × … × m n × ( m 1 × m 3 × … × m n ) − 1 C_2 = r_2 \times m_1 \times m_3 \times \ldots \times m_n \times (m_1 \times m_3 \times \ldots \times m_n)^{-1} C2=r2×m1×m3×…×mn×(m1×m3×…×mn)−1,
- … \ldots …
- C n = r n × m 1 × m 2 × … × m n − 1 × ( m 1 × m 2 × … × m n − 1 ) − 1 C_n = r_n \times m_1 \times m_2 \times \ldots \times m_{n-1} \times (m_1 \times m_2 \times \ldots \times m_{n-1})^{-1} Cn=rn×m1×m2×…×mn−1×(m1×m2×…×mn−1)−1。
当 m 1 , m 2 , … , m n m_1, m_2, \ldots, m_n m1,m2,…,mn 两两互质时,一定存在这样的 C 1 , C 2 , … , C n C_1, C_2, \ldots, C_n C1,C2,…,Cn,因为对应的逆元是一定存在的。且通解为: X = k × lcm + x X = k \times \text{lcm} + x X=k×lcm+x,其中 lcm \text{lcm} lcm因为 m i m_i mi全是质数,所以实际相当于整体的乘积。
总结:在算法竞赛中,中国剩余定理应用的算法流程:
-
计算所有模数的乘积 M = m 1 × m 2 × … × m n M = m_1 \times m_2 \times \ldots \times m_n M=m1×m2×…×mn;
-
计算第 i i i 个方程的 c i = M m i c_i = \frac{M}{m_i} ci=miM;
-
计算 c i c_i ci 在模 m i m_i mi 意义下的逆元 c i − 1 c_i^{-1} ci−1(扩欧算法);
-
最终解为 x = ∑ i = 1 n r i c i c i − 1 ( mod M ) x = \sum_{i=1}^n r_i c_i c_i^{-1} (\text{mod } M) x=∑i=1nricici−1(mod M)。
中国剩余定理OJ示例
模板题曹冲养猪 - 洛谷
P1495 【模板】中国剩余定理(CRT)/ 曹冲养猪 - 洛谷
1634:【例 4】曹冲养猪
CRT的模板题,因为数据量大,需要使用快速幂改装的龟速乘算法来进行最后的求乘。
龟速乘是能一步出结果的乘法,拆分成若干步乘完,这样做能最大程度保证不溢出,但牺牲了速度。
#include<bits/stdc++.h>
using namespace std;typedef long long LL;//扩展欧几里得算法
void exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1; y = 0;return;}LL x1, y1;exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;
}//龟速乘,防溢出
LL qmul(LL a, LL b, LL m) {LL ans = 0;while (b) {if (b % 2)ans = (ans%m + a%m) % m;a = (a + a) % m;b /= 2;}return ans;
}void ac() {LL n;cin >> n;vector<LL>r(n + 1, 0), m(n + 1, 1);for (LL i = 1; i <= n; i++)cin >> m[i] >> r[i];//中国剩余定理LL M = 1;for (auto& x : m)M *= x;LL ans = 0;for (LL i = 1; i <= n; i++) {LL ci = M / m[i];LL x, y;exgcd(ci, m[i], x, y);x = (x % m[i] + m[i]) % m[i];//快速幂改装的鬼速乘算法ans = (ans % M + qmul(qmul(r[i], ci, M), x, M)) % M;}cout << ans << endl;
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;//cin >> T;while (T--)ac();return 0;
}
模板题猜数字 - 洛谷
[P3868 TJOI2009] 猜数字 - 洛谷
中国剩余定理的推广模板。题目要求的是 b i ∣ ( n − a i ) b_i\mid (n-a_i) bi∣(n−ai),即 ( n − a i ) m o d b i = 0 (n-a_i)\bmod b_i=0 (n−ai)modbi=0,将等式拆分开就是 n m o d b i − a i m o d b i = 0 n\bmod b_i-a_i\bmod b_i=0 nmodbi−aimodbi=0,即 n ≡ a i ( m o d b i ) n\equiv a_i\pmod {b_i} n≡ai(modbi)。
因此题目要求的就是同余方程组 n ≡ a i ( m o d b i ) n\equiv a_i\pmod {b_i} n≡ai(modbi)的解, n n n是未知数。
#include<bits/stdc++.h>
using namespace std;typedef long long LL;void exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1; y = 0;return;}LL x1, y1;exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;
}LL qmul(LL a, LL b, LL m) {LL ans = 0;while (b) {if (b % 2)ans = (ans + a) % m;a = (a * 2) % m;b /= 2;}return ans;
}void ac() {int n;cin >> n;vector<LL>a(n + 1, 1), b(a);for (int i = 1; i <= n; i++)cin >> a[i];LL M = 1;for (int i = 1; i <= n; i++)cin >> b[i], M *= b[i];//CRT推广LL ans = 0;for (int i = 1; i <= n; i++) {LL ci = M / b[i];LL x, y;exgcd(ci, b[i], x, y);x = (x % b[i] + b[i]) % b[i];//a[i]为负数,需要放第1个形参,否则会死循环ans = (ans + qmul(qmul(a[i], ci, M), x, M) + M) % M;}cout << ans << endl;
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;//cin >> T;while (T--)ac();return 0;
}
扩展中国剩余定理
中国剩余定理只能处理所有模数两两互质的情况,如果模数不互质,CRT就不适用了。
【扩展中国剩余定理 - EXCRT】
同样是求线性同余方程组:
{ x ≡ r 1 ( mod m 1 ) x ≡ r 2 ( mod m 2 ) ⋮ x ≡ r n ( mod m n ) \begin{cases} x \equiv r_1 (\text{mod } m_1) \\ x \equiv r_2 (\text{mod } m_2) \\ \vdots \\ x \equiv r_n (\text{mod } m_n) \end{cases} ⎩ ⎨ ⎧x≡r1(mod m1)x≡r2(mod m2)⋮x≡rn(mod mn)
思路是将任意整数 X X X 不断迭代成最终解。
这里补充一个方程: x ≡ r 0 ( mod m 0 ) x \equiv r_0 (\text{mod } m_0) x≡r0(mod m0),其中 r 0 = 0 , m 0 = 1 r_0 = 0, m_0 = 1 r0=0,m0=1。
设最终的通解为: r e t = x × m 0 + r 0 ret = x \times m_0 + r_0 ret=x×m0+r0,初始时 m 0 = 1 , r 0 = 0 m_0 = 1, r_0 = 0 m0=1,r0=0。
然后依次迭代后面的方程,使的最终解 r e t ret ret 依次满足后续的方程。
设当前取的方程为: x ≡ r i ( mod m i ) x \equiv r_i (\text{mod } m_i) x≡ri(mod mi),即 r e t = y × m i + r i ret = y \times m_i + r_i ret=y×mi+ri
-
根据两式右边相等: r e t = x × m 0 + r 0 = y × m i + r i ret = x \times m_0 + r_0 = y \times m_i + r_i ret=x×m0+r0=y×mi+ri,
-
于是: x × m 0 − y × m i = r i − r 0 x \times m_0 - y \times m_i = r_i - r_0 x×m0−y×mi=ri−r0,
-
形如 a x + b y = c ax + by = c ax+by=c,其中 a = m 0 , b = m i , c = r i − r 0 a = m_0, b = m_i, c = r_i - r_0 a=m0,b=mi,c=ri−r0
-
根据扩欧算法,可以解出特解 x 0 , y 0 x_0, y_0 x0,y0。(如果无解,算法结束)
-
通解为 x = x 0 + b d × k x = x_0 + \frac{b}{d} \times k x=x0+db×k,代入最终解中: r e t = ( x 0 + b d × k ) × m 0 + r 0 ret = (x_0 + \frac{b}{d} \times k) \times m_0 + r_0 ret=(x0+db×k)×m0+r0
-
整理得: r e t = k ( b d × m 0 ) + ( x 0 × m 0 + r 0 ) ret = k(\frac{b}{d} \times m_0) + (x_0 \times m_0 + r_0) ret=k(db×m0)+(x0×m0+r0),既满足之前的方程,也满足当前方程。
其中新的 m ′ = b d × m 0 , r ′ = x 0 × m 0 + r m' = \frac{b}{d} \times m_0, r' = x_0 \times m_0 + r m′=db×m0,r′=x0×m0+r
- 利用 r e t ret ret的表达式继续迭代下一个方程。直到迭代完所有的方程或某个方程无解为止。
例如:
{ x ≡ 2 ( mod 3 ) x ≡ 3 ( mod 5 ) x ≡ 2 ( mod 7 ) \begin{cases} x \equiv 2(\text{mod } 3) \\ x \equiv 3(\text{mod } 5) \\ x \equiv 2(\text{mod } 7) \end{cases} ⎩ ⎨ ⎧x≡2(mod 3)x≡3(mod 5)x≡2(mod 7)
补上一个方程 x ≡ 0 ( mod 1 ) x \equiv 0 (\text{mod } 1) x≡0(mod 1),最终解 r e t = x × 1 + 0 ret = x \times 1 + 0 ret=x×1+0
迭代第一个方程: r e t = y × 3 + 2 ret = y \times 3 + 2 ret=y×3+2
-
r e t = x × 1 + 0 = y × 3 + 2 ret = x \times 1 + 0 = y \times 3 + 2 ret=x×1+0=y×3+2,即 x − 3 y = 2 x - 3y = 2 x−3y=2,也可看做 x + 3 y = 2 x + 3y = 2 x+3y=2。
-
利用扩欧算法得
{ x = 2 + 3 k y = 0 − k \begin{cases}x=2+3k\\y=0-k\end{cases} {x=2+3ky=0−k
- 此时 r e t = k × 3 + 2 ret=k\times 3+2 ret=k×3+2 满足第一个方程。 k k k是整数,在后续迭代中也可看成新的未知数。
迭代第二个方程: r e t = y × 5 + 3 ret = y \times 5 + 3 ret=y×5+3
-
r e t = x × 3 + 2 = y × 5 + 3 ret = x \times 3 + 2 = y \times 5 + 3 ret=x×3+2=y×5+3,即 3 x + 5 y = 1 3x + 5y = 1 3x+5y=1
-
利用扩欧算法解得:
{ x = 2 + 5 × k y = − 1 − 3 × k \begin{cases} x = 2 + 5 \times k \\ y = -1 - 3 \times k \end{cases} {x=2+5×ky=−1−3×k -
于是将 x x x的通解代入 r e t ret ret 中得: r e t = k × 15 + 8 ret = k \times 15 + 8 ret=k×15+8。此时 r e t ret ret 满足前两个方程。
迭代第三个方程: r e t = y × 7 + 2 ret = y \times 7 + 2 ret=y×7+2
-
r e t = x × 15 + 8 = y × 7 + 2 ret = x \times 15 + 8 = y \times 7 + 2 ret=x×15+8=y×7+2,即 15 x + 7 y = − 6 15x + 7y = -6 15x+7y=−6
-
利用扩欧算法解得:
{ x = − 6 + 7 × k y = 12 − 15 × k \begin{cases} x = -6 + 7 \times k \\ y = 12 - 15 \times k \end{cases} {x=−6+7×ky=12−15×k -
代入 r e t ret ret 中得: r e t = k × 105 − 82 = k × 105 + 23 ret = k \times 105 - 82 = k \times 105 + 23 ret=k×105−82=k×105+23。
此时 r e t ret ret 满足所有方程,并且 23 23 23 为最小非负整数解。
扩展中国剩余定理OJ示例
之前的OJ例如P1495 【模板】中国剩余定理(CRT)/ 曹冲养猪 - 洛谷、[P3868 TJOI2009] 猜数字 - 洛谷也能用扩展中国剩余定理解决。
模板题扩展中国剩余定理(EXCRT) - 洛谷
P4777 【模板】扩展中国剩余定理(EXCRT) - 洛谷
需要注意的点是,构造的方程 a x + b y = c ax+by=c ax+by=c中 c c c需要为最小整数,以及利用 a x + b y = gcd ( a , b ) ax+by=\text{gcd}(a,b) ax+by=gcd(a,b)求方程 a x + b y = c ax+by=c ax+by=c的特解时,需要使用龟速乘算法,这些都是为防止溢出需要做的操作。
#include<bits/stdc++.h>
using namespace std;typedef long long LL;LL exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1; y = 0;return a;}LL x1, y1, d = exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;return d;
}LL qmul(LL a, LL b, LL m) {LL ans = 0;while (b) {if (b & 1)ans = (ans + a) % m;a = (a * 2) % m;b >>= 1;}return ans;
}void ac() {int n;cin >> n;vector<LL>a(n + 1, 0), b(a);for (int i = 1; i <= n; i++)cin >> a[i] >> b[i];//扩展中国剩余定理//ret = k*m+rLL m = 1, r = 0;for (int i = 1; i <= n; i++) {//解不定方程mx+a[i]y=b[i]-rLL x, y, d = exgcd(m, a[i], x, y), c = b[i] - r;c = (c % a[i] + a[i]) % a[i];if (c % d) {cout << -1 << endl;return;}y = a[i] / d;//y是不定方程的另一个解,不重要,于是代码复用//加龟速乘,擦边通关OJx = qmul(x, c / d, y);x = (x % y + y) % y;r = x * m + r;m = y * m;r = (r % m + m) % m;//ret = k*m+r,可以用m对r进行范围限制}cout << r << endl;
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;//cin >> T;while (T--)ac();return 0;
}
NOI2018 屠龙勇士
[P4774 NOI2018] 屠龙勇士 - 洛谷
这种个人感觉是真正的压轴题或准压轴题级别的题,题目又臭又长,还得逐句分析,其中分析还得分阶段,严重的还有分类谈论,每个阶段的分析还得有联系,不然无法流畅地翻译成自然语言或代码。
在实现时往往伴随着知识点的嵌套使用,极其考验综合能力。
读题:
首先是游戏规则:
-
1条龙,打它会使Hp会变负数,然后回血,Hp恢复到0时龙死去,否则不能死去。
例如龙的初始Hp为 5,将血量斩杀到 -4,龙的回血 p = 2 p=2 p=2,回 2 次后龙会死去;或刚好砍到Hp位0,则龙死去。
但如果某几次攻击将血量斩杀到-3,龙回血2次后血量变成1,龙就还能蹦迪。
所以选择的武器需要将龙的血量控制在 p p p的整数倍才能击杀。
-
砍龙需要选剑,砍死龙得新剑。
然后这人想写个脚本之类的东西挂机刷游戏,脚本原理:
-
选择当前有的、伤害不大于龙Hp、同时伤害是最大的剑。没有就选伤害最小的。
例如Hp = 5,剑 a [ i ] = { 1 , 2 , 3 , 4 , 4 , 5 , 6 , 7 , 8 , 9 } a[i]=\{1,2,3,4,4,5,6,7,8,9\} a[i]={1,2,3,4,4,5,6,7,8,9},此时选
a[5]=5
的剑。剑 a [ i ] = 6 , 7 a[i]={6,7} a[i]=6,7,此时选
a[0]=6
的剑。 -
对每条龙,设置固定的攻击次数 x x x。
之后就是解题思路:
-
确定哪条龙用哪把剑。
既然要选剑,则剑应该有序,想到排序。
对每条龙,考虑Hp和奖励,在剑的集合中选不大于龙的血量且伤害最大的,首先应该想到二分优化的枚举,二分的话要考虑数组,数组需要看数据量是否支持。
因为有消耗和奖励,需要频繁插入和删除还有取出指定位置的剑,想到集合
multiset
,插入和取出都是 O ( log n ) O(\log n) O(logn)。可以选择用upper_bound
快速查找,但要注意这个成员函数找的是大于Hp的最小值,所以需要使返回的迭代器减1,边界情况是所有剑都大于Hp。 -
如何确定攻击次数 x x x。
根据题意解读的每条龙的信息
{初始血量Hp,恢复rv,剑的伤害swd,打的次数x,恢复次数y}
,根据游戏规则构造等式 H p − s w d × x + r v × y = 0 Hp-swd\times x+rv\times y=0 Hp−swd×x+rv×y=0,即 s w d × x = r v × y + H p swd\times x=rv\times y+Hp swd×x=rv×y+Hp,看到这个等式想到
s w d × x ÷ r v = y … H p swd\times x\div rv=y\dots Hp swd×x÷rv=y…Hp,2个等式 s w d × x = r v × y + H p swd\times x=rv\times y+Hp swd×x=rv×y+Hp和 s w d × x ÷ r v = y … H p swd\times x\div rv=y\dots Hp swd×x÷rv=y…Hp中, y y y也是未知数,所以根据这2个等式可得出同余方程 s w d × x ≡ H p ( m o d r v ) swd\times x\equiv Hp\pmod {rv} swd×x≡Hp(modrv)。对所有龙的同余方程组:
{ s w d 1 x ≡ H p 1 ( m o d r v 1 ) s w d 2 x ≡ H p 2 ( m o d r v 2 ) … s w d n x ≡ H p n ( m o d r v n ) \begin{cases}swd_1x\equiv Hp_1\pmod {rv_1}\\swd_2x\equiv Hp_2\pmod {rv_2}\\\dots\\swd_nx\equiv Hp_n\pmod {rv_n}\end{cases} ⎩ ⎨ ⎧swd1x≡Hp1(modrv1)swd2x≡Hp2(modrv2)…swdnx≡Hpn(modrvn)
解出 x x x的,即为攻击次数。和之前的同余方程组不同,这个方程组的未知数带系数,使用扩展中国剩余定理时需要注意细节。
-
扩展中国剩余定理解同余方程组:
约定攻击次数 r e t ret ret,恢复力 m = r v m=rv m=rv, r = H p r=Hp r=Hp。
设方程 r e t = k m + r ret=km+r ret=km+r, m = 1 m=1 m=1, r = 0 r=0 r=0。
对任一方程 s w d × x ≡ H p ( m o d r v ) swd\times x\equiv Hp\pmod {rv} swd×x≡Hp(modrv),其实就是 s w d × r e t = y × r v + H p swd\times ret=y\times rv+Hp swd×ret=y×rv+Hp,次数 x x x换成了 r e t ret ret,都表示攻击次数。
联立:
{ r e t = k m + r s w d × r e t = y × r v + H p \begin{cases}ret=km+r\\swd\times ret=y\times rv+Hp\end{cases} {ret=km+rswd×ret=y×rv+Hp
第1个方程左右同时乘 s w d swd swd消去 r e t ret ret,并移项得
s w d × m × k + y × r v = H p − s w d × r swd\times m\times k+y\times rv=Hp-swd\times r swd×m×k+y×rv=Hp−swd×r,此时得到一个不定方程,
对比 A x + B y = C Ax+By=C Ax+By=C, A = s w d × m A=swd\times m A=swd×m, B = r v B=rv B=rv, C = H p − s w d × r C=Hp-swd\times r C=Hp−swd×r。
之后就是熟悉的流程:扩展欧几里得算法求 k k k(或者说 x x x)的通解,代入 r e t = k m + r ret=km+r ret=km+r得到新的 m m m和 r r r,并用它们继续迭代其他方程。
对比之前的扩展中国剩余定理,多个常数系数时只需要联立方程组即可解决。
-
细节问题:
-
防溢出。
在乘法阶段用龟速乘算法代替。同时对方程 s w d × x ≡ H p ( m o d r v ) swd\times x\equiv Hp\pmod {rv} swd×x≡Hp(modrv), s w d × x swd\times x swd×x和 H p Hp Hp可以都同时先模 m m m限制大小,即 ( s w d × x m o d r v ) ≡ ( H p m o d r v ) ( m o d r v ) (swd\times x\bmod rv)\equiv {(Hp\bmod rv)}\pmod {rv} (swd×xmodrv)≡(Hpmodrv)(modrv)。其他地方能取模尽量取模。
-
最终结果 r e t = k m + r ret=km+r ret=km+r中的 r r r可能不是正确答案,而是只是方程组的解,此时还需要 y > 0 y>0 y>0也就是恢复次数也要大于0。
此时对于每条龙的血量 H p Hp Hp和剑的伤害 s w d swd swd,在求攻击次数时需为 ⌈ H p s w d ⌉ \lceil\frac{Hp}{swd}\rceil ⌈swdHp⌉。例如Hp为5的龙,伤害为2的剑,攻击次数需要为 ⌈ 5 2 ⌉ = 3 \lceil\frac{5}{2}\rceil=3 ⌈25⌉=3,才能压血线到负数。
同时需要将每条龙斩到负数Hp的攻击次数并不完全相同,此时需要取最大的次数,取方程中攻击次数大于最大攻击次数的特解,才能保证将所有的龙的血量砍到负数。这里假设所有龙中需要被砍的最多的次数为 max_tm \text{max\_tm} max_tm, max_tm = ⌈ Hp swd ⌉ \text{max\_tm}=\lceil\frac{\text{Hp}}{\text{swd}}\rceil max_tm=⌈swdHp⌉。
对表示攻击次数的方程 r e t = k m + r ret=km+r ret=km+r,需有 r e t ≥ max_tm ret\geq \text{max\_tm} ret≥max_tm。
-
r ≥ max_tm r\geq \text{max\_tm} r≥max_tm,此时 r e t = r ret=r ret=r。
-
r < max_tm r<\text{max\_tm} r<max_tm,此时 r e t = ⌈ max_tm − r m ⌉ × m + r ret=\lceil\frac{\text{max\_tm}-r}{m}\rceil\times m+r ret=⌈mmax_tm−r⌉×m+r。
这里 ⌈ max_tm − r m ⌉ \lceil\frac{\text{max\_tm}-r}{m}\rceil ⌈mmax_tm−r⌉表示需要在最小正整数解的基础上增加的模数 m m m的数量,这样做是为了求能以最低的固定攻击次数将所有的龙都压到负数Hp的特解。
-
-
-
参考程序(变量名和数组名与分析对应):
#include<bits/stdc++.h>
using namespace std;typedef long long LL;LL exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1; y = 0;return a;}LL x1, y1, d = exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;return d;
}LL qmul(LL a, LL b, LL m) {LL ans = 0;while (b) {if (b % 2)ans = (ans + a) % m;a = (a * 2) % m;b /= 2;}return ans;
}void ac() {int N, M;cin >> N >> M;//Hp血量,rv恢复,rd reword奖励,swd sword砍龙时用的剑vector<LL>Hp(N + 1, 0), rv(Hp), rd(Hp),swd(Hp);multiset<LL>st;LL max_tm = 0;//max times,最难砍的龙需要的次数for (int i = 1; i <= N; i++)cin >> Hp[i];for (int i = 1; i <= N; i++)cin >> rv[i];for (int i = 1; i <= N; i++)cin >> rd[i];for (int i = 1; i <= M; i++) {LL x; cin >> x;st.insert(x);}//选择剑for (int i = 1; i <= N; i++) {//二分auto it = st.upper_bound(Hp[i]);if (it != st.begin())--it;swd[i] = *it;max_tm = max(max_tm, (Hp[i] + swd[i] - 1) / swd[i]);st.erase(it);st.insert(rd[i]);}//确定攻击次数//扩展中国剩余定理,初始方程ret=mx+r,m是rv,r是hp,res是攻击次数//对应同余方程(组)swd*x % rv = HpLL m = 1, r = 0;for (int i = 1; i <= N; i++) {//Hp-swd*x+rv*y=0 得到同余方程//联立ret=mx+r和Hp-swd*x+rv*y=0//得到不定方程 (swd-m) x + rv* y= Hp-swd*rLL A = qmul(swd[i] ,m,rv[i]), B = rv[i], C = Hp[i] - swd[i] * r;C = (C % B + B) % B;LL x, y, d = exgcd(A, B, x, y);if (C % d) {cout << -1 << endl;return;}LL g1 = B / d;x = qmul(x,C / d,g1);x = (x % g1 + g1) % g1;r = r + x * m;m = g1 * m;r = (r % m + m) % m;}//最小正整数解可能无法将所有的龙都砍到负数Hp需要求特解if (r < max_tm)r = r + (max_tm - r + m - 1) / m * m;cout << r<<endl;
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;cin >> T;while (T--)ac();return 0;
}