BSGS算法（大步小步法）及其扩展原理及应用
BSGS算法（Baby-step giant-s …
作者： KSkun
中国剩余定理及其扩展原理及应用
中国剩余定理（Chinese remainder theorem）
内容
对于以下一元线性同余方程组：
(S): \begin{cases} x \equiv a_1 \pmod{m_1} \\ x \equiv a_2 \pmod{m_2} \\ \vdots \\ x \equiv a_n \pmod{m_n} \end{cases}
其中整数m_1, m_2, \ldots, m_n两两互质，则对于任意的整数a_1, a_2, \ldots, a_n，方程组 (S) 有解，并且可以用如下方式构造解：
- 设M=\prod_{i=1}^n m_i，M_i = \frac{M}{m_i}
- 设M_i^{-1}为M_i模m_i意义下的逆元，即M_iM_i^{-1} \equiv 1 \pmod{m_i}
- 该方程组的通解为x = kM + \sum_{i=1}^n a_iM_iM_i^{-1} \ (k \in \mathbb{Z})。在模M意义下，该方程组的唯一解为x = \sum_{i=1}^n a_iM_iM_i^{-1}。
例题：中国剩余定理 问题 – 51Nod
按照上述方式构造解即可。代码如下。
// Code by KSkun, 2018/4
#include <cstdio>
typedef long long LL;
inline char fgc() {
static char buf[100000], *p1 = buf, *p2 = buf;
return p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 100000, stdin), p1 == p2) ? EOF : *p1++;
}
inline LL readint() {
register LL res = 0, neg = 1;
char c = fgc();
while(c < '0' || c > '9') {
if(c == '-') neg = -1;
c = fgc();
}
while(c >= '0' && c <= '9') {
res = res * 10 + c - '0';
c = fgc();
}
return res * neg;
}
const int MAXN = 20;
int n, p[MAXN], m[MAXN];
inline LL fpow(LL n, LL k, LL p) {
LL res = 1;
while(k) {
if(k & 1) res = res * n % p;
n = n * n % p;
k >>= 1;
}
return res;
}
inline LL crt() {
LL pall = 1, res = 0;
for(int i = 1; i <= n; i++) {
pall *= p[i];
}
for(int i = 1; i <= n; i++) {
res = (res + m[i] * pall / p[i] % pall
* fpow(pall / p[i], p[i] - 2, p[i]) % pall) % pall;
}
return res;
}
int main() {
n = readint();
for(int i = 1; i <= n; i++) {
p[i] = readint(); m[i] = readint();
}
printf("%lld", crt());
return 0;
}
扩展中国剩余定理（exCRT）
内容
对于CRT的扩展适用于CRT模数非互质的情况。
我们来研究两个同余方程的情况：
\begin{cases} x \equiv b_1 \pmod{a_1} \\ x \equiv b_2 \pmod{a_2} \end{cases} \Rightarrow \begin{cases} x = a_1x_1 + b_1 \\ x = a_2x_2 + b_2 \end{cases}
我们可以联立后面的展开式，得到这个式子a_1x_1 + a_2x_2 = b_2 - b_1。这个式子有解的条件是\mathrm{gcd}(a_1, a_2) | (b_2 - b_1)，如果有解，我们可以利用扩展欧几里得算法（欧几里得算法和扩展欧几里得算法 | KSkun’s Blog）求出x_1的一个解（实际上求出的是方程a_1x_1 + a_2x_2 = \mathrm{gcd}(a_1, a_2)的解，需要转换，即乘一个\frac{b_2-b_1}{\mathrm{gcd}(a_1, a_2)}），回代可以得到x的一个解x'。
这时，我们可以用这个方程x \equiv x' \pmod{\mathrm{lcm}(a_1, a_2)}代替原来的两个方程。
我们可以每次合并两个方程，最终得到一个同余方程，就可以求解了。
例题：2891 — Strange Way to Express Integers
直接利用上述方法解决即可。代码如下。
// Code by KSkun, 2018/4
#include <cstdio>
typedef long long LL;
const int MAXN = 100005;
LL k, a[MAXN], r[MAXN];
inline LL exgcd(LL a, LL b, LL &x, LL &y) {
if(!b) {
x = 1; y = 0; return a;
}
LL res = exgcd(b, a % b, x, y);
LL t = x; x = y; y = t - a / b * y;
return res;
}
inline LL excrt() {
LL A = a[1], R = r[1], x, y;
for(int i = 2; i <= k; i++) {
LL g = exgcd(A, a[i], x, y);
if((r[i] - R) % g) return -1;
x = (r[i] - R) / g * x; x = (x % (a[i] / g) + a[i] / g) % (a[i] / g);
R = A * x + R;
A = A / g * a[i]; R %= A;
}
return R;
}
int main() {
while(scanf("%lld", &k) != EOF) {
for(int i = 1; i <= k; i++) {
scanf("%lld%lld", &a[i], &r[i]);
}
printf("%lld\n", excrt());
}
return 0;
}
参考资料
- 中国剩余定理 – 维基百科，自由的百科全书
- 51nod1079(中国剩余定理) – geloutingyu – 博客园
- 中国剩余定理学习笔记 – MashiroSky – 博客园
- 中国剩余定理与扩展 – CSDN博客
- 《算法竞赛入门经典 训练指南》，刘汝佳，陈锋著，清华大学出版社，9787302291077
![[国家集训队]礼物 题解](https://ksmeow.moe/wp-content/uploads/2018/04/180410a-315x400.jpg "180410a")
![[SDOI2015]星际战争 题解](https://ksmeow.moe/wp-content/uploads/2018/04/180406.jpg "180406")
[SDOI2015]星际战争 题解
题目地址：洛谷：【P3324】[SDOI2015]星际战争 – 洛谷、BZOJ：Problem 3993. — [SDOI2015]星际战争
题目描述
3333年，在银河系的某星球上，X军团和Y军团正在激烈地作战。
在战斗的某一阶段，Y军团一共派遣了N个巨型机器人进攻X军团的阵地，其中第i个巨型机器人的装甲值为Ai。当一个巨型机器人的装甲值减少到0或者以下时，这个巨型机器人就被摧毁了。
X军团有M个激光武器，其中第i个激光武器每秒可以削减一个巨型机器人Bi的装甲值。激光武器的攻击是连续的。
这种激光武器非常奇怪，一个激光武器只能攻击一些特定的敌人。Y军团看到自己的巨型机器人被X军团一个一个消灭，他们急需下达更多的指令。
为了这个目标，Y军团需要知道X军团最少需要用多长时间才能将Y军团的所有巨型机器人摧毁。但是他们不会计算这个问题，因此向你求助。
输入输出格式
输入格式：
第一行，两个整数，N、M。第二行，N个整数，A1、A2…AN。第三行，M个整数，B1、B2…BM。接下来的M行，每行N个整数，这些整数均为0或者1。这部分中的第i行的第j个整数为0表示第i个激光武器不可以攻击第j个巨型机器人，为1表示第i个激光武器可以攻击第j个巨型机器人。
输出格式：
一行，一个实数，表示X军团要摧毁Y军团的所有巨型机器人最少需要的时间。输出结果与标准答案的绝对误差不超过10-3即视为正确。
输入输出样例
输入样例#1：
2 2 3 10 4 6 0 1 1 1
输出样例#1：
1.300000
说明
【样例说明1】
战斗开始后的前0.5秒，激光武器1攻击2号巨型机器人，激光武器2攻击1号巨型机器人。1号巨型机器人被完全摧毁，2号巨型机器人还剩余8的装甲值；
接下来的0.8秒，激光武器1、2同时攻击2号巨型机器人。2号巨型机器人被完全摧毁。
对于全部的数据，1<=N, M<=50，1<=Ai<=10510^5105 ，1<=Bi<=1000，输入数据保证X军团一定能摧毁Y军团的所有巨型机器人
题解
这个问题本身很难办，因为一个武器无法同时攻击多个目标，这一点在常规的网络流算法中难以控制。既然如此，我们考虑转换为判定问题，判断所给时间内是否能打败。这个我们可以用最大流跑，建一个源→武器→机器人→汇的图，其中源→武器的容量为武器在这段时间内的攻击输出（即t \cdot B_i），机器人→汇的容量为机器人的装甲值。只要汇点可以满流，就说明可以打败。
至于精度问题，要求1e-3的误差，那么扩大1e4倍即可解决问题。
代码
// Code by KSkun, 2018/4
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <queue>
typedef long long LL;
const int MAXN = 105, MAXM = 500005;
const LL INF = 1e15;
const double EPS = 1e-6;
struct Edge {
int to;
LL cap;
int nxt;
} gra[MAXM << 1];
int head[MAXN], tot;
inline void addedge(int u, int v, LL cap) {
gra[tot] = Edge {v, cap, head[u]}; head[u] = tot++;
gra[tot] = Edge {u, 0, head[v]}; head[v] = tot++;
}
inline void init() {
memset(head, -1, sizeof(head));
tot = 0;
}
int level[MAXN];
bool vis[MAXN];
std::queue<int> que;
inline bool bfs(int s, int t) {
memset(level, -1, sizeof(level));
level[s] = 0; que.push(s);
while(!que.empty()) {
int u = que.front(); que.pop();
for(int i = head[u]; ~i; i = gra[i].nxt) {
int v = gra[i].to;
if(level[v] == -1 && gra[i].cap > 0) {
level[v] = level[u] + 1;
que.push(v);
}
}
}
return level[t] != -1;
}
inline LL dfs(int u, int t, LL left) {
if(left == 0 || u == t) return left;
vis[u] = true;
LL flow = 0;
for(int i = head[u]; ~i; i = gra[i].nxt) {
int v = gra[i].to;
if(!vis[v] && gra[i].cap > 0 && level[v] == level[u] + 1) {
int d = dfs(v, t, std::min(left, gra[i].cap));
if(d > 0) {
gra[i].cap -= d; gra[i ^ 1].cap += d;
flow += d; left -= d;
if(left == 0) {
level[u] = -1;
return flow;
}
}
}
}
return flow;
}
inline LL dinic(int s, int t) {
LL flow = 0;
while(bfs(s, t)) {
memset(vis, 0, sizeof(vis));
LL f;
while(f = dfs(s, t, INF)) {
flow += f;
}
}
return flow;
}
LL n, m, a[MAXN], b[MAXN], S, T, sum;
bool canatt[MAXN][MAXN];
inline bool check(LL mid) {
init();
for(int i = 1; i <= n; i++) {
addedge(S, i, a[i]);
}
for(int i = 1; i <= m; i++) {
addedge(i + n, T, b[i] * mid);
}
for(int i = 1; i <= m; i++) {
for(int j = 1; j <= n; j++) {
if(canatt[i][j]) addedge(j, i + n, INF);
}
}
LL flow = dinic(S, T);
return flow >= sum;
}
int main() {
scanf("%lld%lld", &n, &m);
S = n + m + 1; T = S + 1;
for(int i = 1; i <= n; i++) {
scanf("%lld", &a[i]);
a[i] *= 1e4;
sum += a[i];
}
for(int i = 1; i <= m; i++) {
scanf("%lld", &b[i]);
}
for(int i = 1; i <= m; i++) {
for(int j = 1; j <= n; j++) {
scanf("%d", &canatt[i][j]);
}
}
LL l = 0, r = 1e9, mid;
while(r - l > 1) {
mid = (l + r) >> 1;
if(check(mid)) r = mid; else l = mid;
}
printf("%.4lf", r / 1e4);
return 0;
}
![[USACO15DEC]最大流Max Flow 题解](https://ksmeow.moe/wp-content/uploads/2018/04/180405c-696x400.jpg "180405c")
![[SCOI2007]修车 题解](https://ksmeow.moe/wp-content/uploads/2018/04/180405b-338x400.jpg "180405b")
![[洛谷3380]【模板】二逼平衡树（树套树） 题解](https://ksmeow.moe/wp-content/uploads/2018/04/180405a.jpg "180405a")
[洛谷3380]【模板】二逼平衡树（树套树） 题解
题目地址：洛谷：【P3380】【模板】二逼平衡树（树套树） – 洛谷、BZOJ：Problem 3196. — Tyvj 1730 二逼平衡树
题目描述
您需要写一种数据结构（可参考题目标题），来维护一个有序数列，其中需要提供以下操作：
- 查询k在区间内的排名
- 查询区间内排名为k的值
- 修改某一位值上的数值
- 查询k在区间内的前驱(前驱定义为严格小于x，且最大的数，若不存在输出-2147483647)
- 查询k在区间内的后继(后继定义为严格大于x，且最小的数，若不存在输出2147483647)
输入输出格式
输入格式：
第一行两个数 n,m 表示长度为n的有序序列和m个操作
第二行有n个数，表示有序序列
下面有m行，opt表示操作标号
若opt=1 则为操作1，之后有三个数l,r,k 表示查询k在区间[l,r]的排名
若opt=2 则为操作2，之后有三个数l,r,k 表示查询区间[l,r]内排名为k的数
若opt=3 则为操作3，之后有两个数pos,k 表示将pos位置的数修改为k
若opt=4 则为操作4，之后有三个数l,r,k 表示查询区间[l,r]内k的前驱
若opt=5 则为操作5，之后有三个数l,r,k 表示查询区间[l,r]内k的后继
输出格式：
对于操作1,2,4,5各输出一行，表示查询结果
输入输出样例
输入样例#1：
9 6 4 2 2 1 9 4 0 1 1 2 1 4 3 3 4 10 2 1 4 3 1 2 5 9 4 3 9 5 5 2 8 5
输出样例#1：
2 4 3 4 9
说明
时空限制：2s，128M
n,m≤5*10^4保证有序序列所有值在任何时刻满足[0,10^8]
题目来源：bzoj3196 / Tyvj1730 二逼平衡树，在此鸣谢
此数据为洛谷原创。（特别提醒：此数据不保证操作5、6一定存在，故请务必考虑不存在的情况）
题解
我们考虑使用线段树套Splay，即对于线段树中的每一个区间，建一棵Splay维护其中的元素，每次按照线段树的操作依次更新涉及区间。
对于操作1，当区间割裂开的时候不好处理，我们考虑计算严格小于这个数的个数，割裂的两区间的结果加起来就是答案，而排名即这一结果+1。
对于操作2，不好操作，我们二分这个排名为k的数是哪个，因为一旦比这个数字大，肯定大于这个数字的排名，虽然这个数字的排名可能会小于k，但这并不影响二分的正确性。
对于操作3，erase后insert即可。
对于操作4、5，从割裂区间中的结果取min/max即可。
由于自带大常数，只能在O2的帮助下A本题。
代码
// Code by KSkun, 2018/4
#include <cstdio>
#include <cstring>
inline int min(int a, int b) {
return a < b ? a : b;
}
inline int max(int a, int b) {
return a > b ? a : b;
}
inline char fgc() {
static char buf[100000], *p1 = buf, *p2 = buf;
return p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 100000, stdin), p1 == p2) ? EOF : *p1++;
}
inline int readint() {
register int res = 0, neg = 1;
register char c = fgc();
while(c < '0' || c > '9') {
if(c == '-') neg = -1;
c = fgc();
}
while(c >= '0' && c <= '9') {
res = res * 10 + c - '0';
c = fgc();
}
return res * neg;
}
const int MAXN = 50005, INF = 2147483647;
struct Node {
int fa, ch[2], val, siz, cnt;
} tr[MAXN * 30];
int tot = 0, sta[MAXN], stop = 0;
inline int newnode() {
register int p;
if(stop) p = sta[--stop];
else p = ++tot;
memset(tr + p, 0, sizeof(Node));
return p;
}
inline void delnode(int p) {
sta[stop++] = p;
}
struct Splay {
int rt;
Splay() {
rt = 0;
}
inline void update(int p) {
register int lch = tr[p].ch[0], rch = tr[p].ch[1];
tr[p].siz = tr[lch].siz + tr[rch].siz + tr[p].cnt;
}
inline bool isleft(int p) {
return tr[tr[p].fa].ch[0] == p;
}
inline void rotate(int p) {
register bool t = !isleft(p); register int fa = tr[p].fa, ffa = tr[fa].fa;
tr[p].fa = ffa; if(ffa) tr[ffa].ch[!isleft(fa)] = p;
tr[fa].ch[t] = tr[p].ch[!t]; tr[tr[fa].ch[t]].fa = fa;
tr[p].ch[!t] = fa; tr[fa].fa = p;
update(fa);
if(!tr[p].fa) rt = p;
}
inline void splay(int p, int tar) {
for(register int fa = tr[p].fa; fa != tar; rotate(p), fa = tr[p].fa) {
if(tr[fa].fa != tar) rotate(isleft(fa) == isleft(p) ? fa : p);
}
update(p);
}
inline void insert(int v) {
if(!rt) {
rt = newnode();
tr[rt].val = v;
tr[rt].siz = tr[rt].cnt = 1;
return;
}
register int p = rt, fa = 0;
for(;;) {
if(v == tr[p].val) {
tr[p].cnt++;
update(p); update(fa);
splay(p, 0);
return;
}
fa = p;
p = tr[p].ch[tr[p].val < v];
if(!p) {
p = newnode();
tr[p].val = v;
tr[p].siz = tr[p].cnt = 1;
tr[p].fa = fa;
tr[fa].ch[tr[fa].val < v] = p;
update(fa);
splay(p, 0);
return;
}
}
}
inline int queryrk(int v) {
register int p = rt, res = 0;
for(;;) {
if(v < tr[p].val) {
p = tr[p].ch[0];
} else {
res += tr[tr[p].ch[0]].siz;
if(v == tr[p].val) {
splay(p, 0);
return res + 1;
}
res += tr[p].cnt;
p = tr[p].ch[1];
}
}
}
inline int queryn(int rk) {
register int p = rt;
for(;;) {
if(tr[p].ch[0] && rk <= tr[tr[p].ch[0]].siz) {
p = tr[p].ch[0];
} else {
rk -= tr[tr[p].ch[0]].siz;
if(rk <= tr[p].cnt) {
return tr[p].val;
}
rk -= tr[p].cnt;
p = tr[p].ch[1];
}
}
}
inline int querypre() {
register int p = tr[rt].ch[0];
while(tr[p].ch[1]) p = tr[p].ch[1];
return p;
}
inline int querynxt() {
register int p = tr[rt].ch[1];
while(tr[p].ch[0]) p = tr[p].ch[0];
return p;
}
inline void erase(int v) {
queryrk(v);
if(tr[rt].cnt > 1) {
tr[rt].cnt--;
update(rt);
return;
}
if(!tr[rt].ch[0]) {
delnode(rt);
rt = tr[rt].ch[1];
tr[rt].fa = 0;
return;
}
if(!tr[rt].ch[1]) {
delnode(rt);
rt = tr[rt].ch[0];
tr[rt].fa = 0;
return;
}
register int ort = rt, lmx = querypre();
splay(lmx, 0);
tr[rt].ch[1] = tr[ort].ch[1];
tr[tr[rt].ch[1]].fa = rt;
delnode(ort);
update(rt);
}
} spt[MAXN << 2];
int n, m, a[MAXN];
#define lch o << 1
#define rch (o << 1) | 1
#define mid ((l + r) >> 1)
inline void build(int o, int l, int r) {
for(register int i = l; i <= r; i++) {
spt[o].insert(a[i]);
}
if(l == r) return;
build(lch, l, mid);
build(rch, mid + 1, r);
}
inline int querybef(int o, int l, int r, int ll, int rr, int k) {
if(l == ll && r == rr) {
spt[o].insert(k);
register int res = spt[o].queryrk(k) - 1;
spt[o].erase(k);
return res;
}
if(rr <= mid) {
return querybef(lch, l, mid, ll, rr, k);
} else if(ll > mid) {
return querybef(rch, mid + 1, r, ll, rr, k);
} else {
return querybef(lch, l, mid, ll, mid, k) + querybef(rch, mid + 1, r, mid + 1, rr, k);
}
}
inline int queryk(int ll, int rr, int k) {
register int l = 0, r = 100000001;
while(r - l > 1) {
register int rk = querybef(1, 1, n, ll, rr, mid) + 1;
if(rk <= k) l = mid; else r = mid;
}
return l;
}
inline void modify(int o, int l, int r, int x, int w) {
spt[o].erase(a[x]);
spt[o].insert(w);
if(l == r) return;
if(x <= mid) modify(lch, l, mid, x, w);
else modify(rch, mid + 1, r, x, w);
}
inline int querypre(int o, int l, int r, int ll, int rr, int k) {
if(l == ll && r == rr) {
spt[o].insert(k);
spt[o].queryrk(k);
register int res = tr[spt[o].querypre()].val;
spt[o].erase(k);
return res ? res : -INF;
}
if(rr <= mid) {
return querypre(lch, l, mid, ll, rr, k);
} else if(ll > mid) {
return querypre(rch, mid + 1, r, ll, rr, k);
} else {
return max(querypre(lch, l, mid, ll, mid, k), querypre(rch, mid + 1, r, mid + 1, rr, k));
}
}
inline int querynxt(int o, int l, int r, int ll, int rr, int k) {
if(l == ll && r == rr) {
spt[o].insert(k);
spt[o].queryrk(k);
register int res = tr[spt[o].querynxt()].val;
spt[o].erase(k);
return res ? res : INF;
}
if(rr <= mid) {
return querynxt(lch, l, mid, ll, rr, k);
} else if(ll > mid) {
return querynxt(rch, mid + 1, r, ll, rr, k);
} else {
return min(querynxt(lch, l, mid, ll, mid, k), querynxt(rch, mid + 1, r, mid + 1, rr, k));
}
}
int op, l, r, x, k;
int main() {
n = readint(); m = readint();
for(register int i = 1; i <= n; i++) {
a[i] = readint();
}
build(1, 1, n);
while(m--) {
op = readint();
if(op != 3) {
l = readint(); r = readint();
} else {
x = readint();
}
k = readint();
switch(op) {
case 1:
printf("%d\n", querybef(1, 1, n, l, r, k) + 1); break;
case 2:
printf("%d\n", queryk(l, r, k)); break;
case 3:
modify(1, 1, n, x, k); a[x] = k; break;
case 4:
printf("%d\n", querypre(1, 1, n, l, r, k)); break;
case 5:
printf("%d\n", querynxt(1, 1, n, l, r, k));
}
}
return 0;
}
![[JSOI2011]任务调度 题解](https://ksmeow.moe/wp-content/uploads/2018/04/180404c-674x400.jpg "180404c")
