随机化杂记
随机化杂记
模拟退火
主要参数:初温 \(T_0\) ,末温 \(T_k\) ,当前温度 \(T\) ,温度变动量 \(\Delta\) ,当前最优解 \(E_0\) ,新解 \(E\) ,上一个被接受的解 \(E_1\) ,解变动量 \(\Delta E\) 。
过程:
- 从 \(T = T_0\) 开始,每次乘上 \(\Delta\) ,当 \(T < T_k\) 时退出。
- 每次在 \(E_1\) 的基础上扰动产生新解 \(E\) ,其中扰动大小随温度的降低而变小。
- Metropolis 准则:用于判断是否接受这个新解。若 \(E\) 优于 \(E_0\) 则一定接受,否则按 \(\exp(\frac{-\Delta E}{T})\) 的概率接受。
一般取 \(T\) 为 \(5000\) 左右,\(\Delta\) 为一个略 \(< 1\) 的常数,\(T_k\) 一般取 \(10^{-9}\) 左右。
用随机集合覆盖目标元素
CF1305F Kuroni and the Punishment
给定 \(n\) 个数。每次可以选择将一个数 \(+1\) 或 \(-1\) 。求至少多少次操作使得整个序列都是正数且全部元素的 \(\gcd>1\) 。
\(n \leq 2 \times 10^5\) ,\(a_i \leq 10^{12}\)
很容易想到一个上界为奇数的数量,构造就是把所有奇数 \(+1\) ,这样 \(\gcd = 2\) 合法。
注意到:必有一半以上元素至多操作一次,证明考虑反证法,若不成立则操作次数超过 \(n\) 显然不优。于是可以类似绝对众数,考虑随机一个数对它进行 \(+1\) 、\(-1\) 或不操作。此时就可以固定一个数,显然钦定 \(\gcd\) 为它的质因子是最优的,于是枚举它的质因子计算操作次数即可。
一次随机到正确的数的概率是 \(\frac{1}{2}\) ,随机约 \(20\) 次即可通过。
#include <bits/stdc++.h>
typedef long long ll;
using namespace std;
const int N = 2e5 + 7, V = 1e6 + 7;
ll a[N];
int pri[V];
bool isp[V];
mt19937 myrand(time(0));
int n, pcnt;
inline void sieve() {
memset(isp, true, sizeof(isp));
isp[1] = false;
for (int i = 2; i < V; ++i) {
if (isp[i])
pri[++pcnt] = i;
for (int j = 1; j <= pcnt && i * pri[j] < V; ++j) {
isp[i * pri[j]] = false;
if (!(i % pri[j]))
break;
}
}
}
inline int calc(ll k) {
ll res = 0;
for (int i = 1; i <= n; ++i) {
if (a[i] < k)
res += k - a[i];
else
res += min(a[i] % k, k - a[i] % k);
if (res > n)
return n;
}
return res;
}
inline int solve(ll k) {
int res = n;
for (int i = 1; 1ll * pri[i] * pri[i] <= k; ++i)
if (!(k % pri[i])) {
res = min(res, calc(pri[i]));
while (!(k % pri[i]))
k /= pri[i];
}
if (k > 1)
res = min(res, calc(k));
return res;
}
signed main() {
sieve();
scanf("%d", &n);
int ans = 0;
for (int i = 1; i <= n; ++i)
scanf("%lld", a + i), ans += a[i] & 1;
for (int i = 1; i <= 20; ++i) {
int id = myrand() % n + 1;
ans = min(ans, min(solve(a[id]), solve(a[id] + 1)));
if (a[id] > 1)
ans = min(ans, solve(a[id] - 1));
}
printf("%d", ans);
return 0;
}
CF1155F Delivery Oligopoly
给出一张边双连通图,去掉最多的边使得剩下的图仍是边双连通图。
\(n \leq 14\) ,\(m \leq \frac{n(n - 1)}{2}\)
考虑一个生成双联通子图的方法:跑一遍 Tarjan,初始时只选择所有树边,每次当一个点所有儿子的 \(low\) 比其 \(dfn\) 要大的时候,就加入一条子树内往上连的尽量高的边。
但是这不能保证是一个最优解,考虑每次随机一个点为根 Tarjan,然后每次打乱边的顺序,多跑几次即可。
#include <bits/stdc++.h>
using namespace std;
const int N = 15;
struct Graph {
vector<int> e[N];
inline void insert(int u, int v) {
e[u].emplace_back(v);
}
} G;
vector<pair<int, int> > ans, vec;
pair<int, int> up[N];
int dfn[N], low[N], low2[N];
mt19937 myrand(time(0));
int n, m, dfstime;
void Tarjan(int u, int f) {
shuffle(G.e[u].begin(), G.e[u].end(), myrand);
dfn[u] = low[u] = low2[u] = ++dfstime;
for (int v : G.e[u]) {
if (!dfn[v]) {
vec.emplace_back(u, v), Tarjan(v, u);
if (low[u] > low[v])
up[u] = up[v];
low[u] = min(low[u], low[v]);
low2[u] = min(low2[u], low2[v]);
} else if (v != f) {
if (dfn[v] < low[u])
up[u] = make_pair(u, v);
low[u] = min(low[u], dfn[v]);
}
}
if (dfn[u] == low2[u] && f)
vec.emplace_back(up[u]), low2[u] = low[u];
}
signed main() {
scanf("%d%d", &n, &m);
for (int i = 1; i <= m; ++i) {
int u, v;
scanf("%d%d", &u, &v);
G.insert(u, v), G.insert(v, u);
ans.emplace_back(u, v);
}
for (int i = 1; i <= 1e3; ++i) {
memset(dfn + 1, 0, sizeof(int) * n);
memset(low + 1, 0, sizeof(int) * n);
memset(low2 + 1, 0, sizeof(int) * n);
vec.clear(), dfstime = 0;
Tarjan(myrand() % n + 1, 0);
if (vec.size() < ans.size())
ans = vec;
}
printf("%d\n", ans.size());
for (auto it : ans)
printf("%d %d\n", it.first, it.second);
return 0;
}
P1224 [NOI2013] 向量内积
给出 \(n\) 个 \(d\) 维向量,找到两个向量满足内积为 \(k\) 的倍数,或报告无解。
\(k \in \{ 2, 3 \}\)
- Task 1:\(n \leq 2 \times 10^4\) ,\(d \leq 100\)
- Taks 2:\(n \leq 10^5\) ,\(d \leq 30\)
先考虑 \(k = 2\) 的情况,若不存在合法的 \((i, j)\) ,则必然有 \(\sum_k a_{i, k} \times a_{j, k} = 1\) 。
考虑随机打乱所有向量,对每个 \(i\) 计算第 \(i\) 个向量与前 \(i - 1\) 个的内积和。由于内积在 \(\bmod 2\) 意义下取值只有 \(0\) 和 \(1\) ,考虑判定这个内积和模 \(k\) 是否为 \((i - 1) \bmod k\) ,若不等于则说明 \(i\) 一定与前面某个向量内积为 \(k\) 的倍数,单次时间复杂度 \(O(nd)\) 。每次正确率至少为 \(\frac{1}{2}\) ,随机多次即可。
接下来考虑 \(k = 3\) 的情况,考虑计算内积的平方。在 \(\bmod 3\) 意义下,若内积非 \(0\) ,则内积的平方必然为 \(1\) 。
由于 \((\sum_{k = 1}^d a_{i, k} \times a_{j, k})^2 = \sum_{x, y} a_{i, x} \times a_{i, y} \times a_{j, x} \times a_{j, y}\) ,因此这可以转化为长度为 \(d^2\) 的向量的内积,单次时间复杂度 \(O(n d^2)\) 。
为了实现方便,其实 \(k = 2\) 的时候也可以用 \(O(n d^2)\) 的算法判定。
#include <bits/stdc++.h>
using namespace std;
const int N = 1e5 + 7, D = 1e2 + 7;
int a[N][D], s[D][D];
mt19937 myrand(time(0));
int n, d, k;
inline int query(int x) {
int res = 0;
for (int i = 1; i <= d; ++i)
for (int j = 1; j <= d; ++j)
res = (res + a[x][i] * a[x][j] * s[i][j]) % k, s[i][j] = (s[i][j] + a[x][i] * a[x][j]) % k;
return res;
}
inline bool check(int x, int y) {
int res = 0;
for (int i = 1; i <= d; ++i)
res = (res + a[x][i] * a[y][i]) % k;
return !(res % k);
}
signed main() {
scanf("%d%d%d", &n, &d, &k);
for (int i = 1; i <= n; ++i)
for (int j = 1; j <= d; ++j)
scanf("%d", a[i] + j), a[i][j] %= k;
vector<int> id(n);
iota(id.begin(), id.end(), 1);
while ((double) clock() / CLOCKS_PER_SEC < 0.8) {
shuffle(id.begin(), id.end(), myrand);
memset(s, 0, sizeof(s));
for (int i = 0; i < n; ++i) {
int res = query(id[i]);
if (res == i % k)
continue;
for (int j = 0; j < i; ++j)
if (check(id[i], id[j]))
return printf("%d %d\n", min(id[i], id[j]), max(id[i], id[j])), 0;
}
}
puts("-1 -1");
return 0;
}
随机哈希映射
P5227 [AHOI2013] 连通图
给定一张无向图,每次询问删除给出的 \(c_i\) 条边后整张图是否连通,询问之间独立。
\(n, q \leq 2 \times 10^5\) ,\(m \leq 2 \times 10^5\) ,\(c_i \leq 4\)
构建 dfs 树,考虑给每条非树边一个随机权值,每条树边的权值为覆盖它的非树边的权值异或和,可以用树上差分解决。
那么若删除的边的一个子集异或和为 \(0\) ,则大概率说明一条树边与覆盖它的非树边都被删除,即原图不连通。
直接上线性基,时间复杂度 \(O(q \log V)\) 。
#include <bits/stdc++.h>
typedef unsigned int uint;
using namespace std;
const int N = 1e5 + 7, M = 2e5 + 7, B = 32;
struct Graph {
vector<pair<int, int> > e[N];
inline void insert(int u, int v, int w) {
e[u].emplace_back(v, w);
}
} G;
struct LinearBasis {
uint p[B];
inline LinearBasis() {
memset(p, 0, sizeof(p));
}
inline bool insert(uint k) {
for (int i = B - 1; ~i; --i)
if (k >> i & 1) {
if (p[i])
k ^= p[i];
else
return p[i] = k, true;
}
return false;
}
};
uint d[N], val[M];
bool vis[N];
mt19937 myrand(time(0));
int n, m, q;
void dfs(int u, int f) {
vis[u] = true;
for (auto it : G.e[u]) {
int v = it.first, id = it.second;
if (!vis[v])
dfs(v, u), d[u] ^= d[v], val[id] = d[v];
else if (v != f && !val[id])
val[id] = myrand(), d[u] ^= val[id], d[v] ^= val[id];
}
}
signed main() {
scanf("%d%d", &n, &m);
for (int i = 1; i <= m; ++i) {
int u, v;
scanf("%d%d", &u, &v);
G.insert(u, v, i), G.insert(v, u, i);
}
dfs(1, 0);
scanf("%d", &q);
while (q--) {
int num;
scanf("%d", &num);
LinearBasis lb;
bool flag = true;
while (num--) {
int x;
scanf("%d", &x);
flag &= lb.insert(val[x]);
}
puts(flag ? "Connected" : "Disconnected");
}
return 0;
}
UOJ207. 共价大爷游长沙
动态维护一棵树和二元组集合,需要支持:
- 删边并加边,保证图仍为树。
- 加入/删除某个二元组。
- 给出一条边 \(e\) ,判定是否对于集合中的每个二元组 \((s, t)\) ,\(e\) 都在 \(s, t\) 的路径上。
\(n \leq 10^5\) ,\(m \leq 3 \times 10^5\)
首先不难想到用 LCT 动态维护树的形态,考虑对每条边用哈希维护经过这条边的二元组集合。
给每个二元组随机赋上一个权值,这样一条边的哈希值就可以定义为子树异或和,加入/删除某个二元组就是单点异或。
时间复杂度 \(O((n + m) \log n)\) 。
#include <bits/stdc++.h>
typedef unsigned long long ull;
using namespace std;
const int N = 1e5 + 7, M = 3e5 + 7;
struct Node {
ull h;
int x, y;
} nd[M];
mt19937_64 myrand(time(0));
int testid, n, m, tot;
namespace LCT {
ull val[N], s[N], s2[N];
int ch[N][2], fa[N], rev[N], sta[N];
inline bool isroot(int x) {
return x != ch[fa[x]][0] && x != ch[fa[x]][1];
}
inline int dir(int x) {
return x == ch[fa[x]][1];
}
inline void pushup(int x) {
s[x] = val[x] ^ s2[x];
if (ch[x][0])
s[x] ^= s[ch[x][0]];
if (ch[x][1])
s[x] ^= s[ch[x][1]];
}
inline void reverse(int x) {
swap(ch[x][0], ch[x][1]), rev[x] ^= 1;
}
inline void pushdown(int x) {
if (rev[x]) {
if (ch[x][0])
reverse(ch[x][0]);
if (ch[x][1])
reverse(ch[x][1]);
rev[x] = 0;
}
}
inline void rotate(int x) {
int y = fa[x], z = fa[y], d = dir(x);
if (!isroot(y))
ch[z][dir(y)] = x;
fa[x] = z, ch[y][d] = ch[x][d ^ 1];
if (ch[x][d ^ 1])
fa[ch[x][d ^ 1]] = y;
ch[x][d ^ 1] = y, fa[y] = x;
pushup(y), pushup(x);
}
inline void splay(int x) {
int top = 0;
sta[++top] = x;
for (int cur = x; !isroot(cur); cur = fa[cur])
sta[++top] = fa[cur];
while (top)
pushdown(sta[top--]);
for (int f = fa[x]; !isroot(x); rotate(x), f = fa[x])
if (!isroot(f))
rotate(dir(f) == dir(x) ? f : x);
}
inline void access(int x) {
for (int y = 0; x; y = x, x = fa[x])
splay(x), s2[x] ^= s[y] ^ s[ch[x][1]], ch[x][1] = y, pushup(x);
}
inline void makeroot(int x) {
access(x), splay(x), reverse(x);
}
inline void split(int x, int y) {
makeroot(x), access(y), splay(y);
}
inline int findroot(int x) {
access(x), splay(x);
while (ch[x][0])
x = ch[x][0];
return x;
}
inline void link(int x, int y) {
if (findroot(x) != findroot(y))
makeroot(x), fa[x] = y, s2[y] ^= s[x], pushup(y);
}
inline void cut(int x, int y) {
if (findroot(x) != findroot(y))
return;
split(x, y);
if (ch[y][0] == x && !ch[x][1])
fa[x] = ch[y][0] = 0, pushup(y);
}
inline void update(int x, ull k) {
makeroot(x), val[x] ^= k, pushup(x);
}
inline ull query(int x, int y) {
return split(x, y), s2[y] ^ val[y];
}
} // namespace LCT
signed main() {
scanf("%d%d%d", &testid, &n, &m);
for (int i = 1; i < n; ++i) {
int u, v;
scanf("%d%d", &u, &v);
LCT::link(u, v);
}
ull all = 0;
while (m--) {
int op;
scanf("%d", &op);
if (op == 1) {
int x, y, u, v;
scanf("%d%d%d%d", &x, &y, &u, &v);
LCT::cut(x, y), LCT::link(u, v);
} else if (op == 2) {
all ^= (nd[++tot].h = myrand());
scanf("%d%d", &nd[tot].x, &nd[tot].y);
LCT::update(nd[tot].x, nd[tot].h), LCT::update(nd[tot].y, nd[tot].h);
} else if (op == 3) {
int x;
scanf("%d", &x);
all ^= nd[x].h, LCT::update(nd[x].x, nd[x].h), LCT::update(nd[x].y, nd[x].h);
} else {
int x, y;
scanf("%d%d", &x, &y);
puts(LCT::query(x, y) == all ? "YES" : "NO");
}
}
return 0;
}
P10785 [NOI2024] 集合
定义:
- 基本集合:大小为 \(3\) ,元素 \(\in [1, m]\) 。
- 集合序列:由基本元素组成的序列。
- 等价:定义两个长度为 \(n\) 的集合序列 \(A, B\) 等价,当且仅当存在一个置换 \(p_{1 \sim m}\) ,满足 \(A_i\) 的元素经过置换后即为 \(B_i\) 。
给出两个长度为 \(n\) 的集合序列 \(A, B\) ,\(q\) 次询问 \(A_{l \sim r}\) 与 \(B_{l \sim r}\) 是否等价。
\(n \leq 2 \times 10^5\) ,\(m \leq 6 \times 10^5\) ,\(q \leq 10^6\)
考虑如何判断两个集合序列是否等价,等价应当要找到一置换 \(p_{1 \sim m}\) ,满足 \(i\) 在 \(A\) 中的出现位置集合与 \(p_i\) 在 \(B\) 中的出现位置集合相等。
考虑随机哈希映射,将每个位置赋上一个随机哈希值,一个数的出现位置集合映射为其出现位置上哈希值的异或和 \(ha_i, hb_i\) ,则判定条件转化为存在一个置换 \(p_{1 \sim m}\) ,满足 \(ha_i = hb_{p_i}\) 。
考虑再次映射,将一个集合序列的哈希值设为 \(1 \sim m\) 的哈希值之和,那么判定条件转化为两个集合哈希值相等。
不难发现当 \(r\) 增加时,\(l\) 是不降的,于是考虑双指针对每个 \(r\) 求出 \(l\) 的范围,不难做到在线 \(O(1)\) 修改某数的哈希值,进而维护集合序列的哈希值。
时间复杂度 \(O(n + q)\) 。
#include <bits/stdc++.h>
typedef unsigned long long ull;
using namespace std;
const int N = 6e5 + 7, S = 3e7 + 7;
ull val[N], hsa[N], hsb[N];
int a[N], b[N], mxl[N];
mt19937_64 myrand(time(0));
int n, m, q;
template <class T = int>
inline T read() {
char c = getchar();
bool sign = (c == '-');
while (c < '0' || c > '9')
c = getchar(), sign |= (c == '-');
T x = 0;
while ('0' <= c && c <= '9')
x = (x << 1) + (x << 3) + (c & 15), c = getchar();
return sign ? ~x + 1 : x;
}
inline ull calc(ull x) {
x = x * x * x;
x ^= x >> 33;
x *= 2398562385683465ull;
x ^= x >> 17;
x = x * (x - 1);
return x;
}
signed main() {
scanf("%d%d%d", &n, &m, &q);
for (int i = 1; i <= n * 3; ++i)
scanf("%d", a + i);
for (int i = 1; i <= n * 3; ++i)
scanf("%d", b + i);
ull suma = 0, sumb = 0;
for (int i = 1, j = 1; i <= n; ++i) {
val[i] = myrand();
for (int k = (i - 1) * 3 + 1; k <= i * 3; ++k) {
suma -= calc(hsa[a[k]]), suma += calc(hsa[a[k]] ^= val[i]);
sumb -= calc(hsb[b[k]]), sumb += calc(hsb[b[k]] ^= val[i]);
}
for (; suma != sumb; ++j)
for (int k = (j - 1) * 3 + 1; k <= j * 3; ++k) {
suma -= calc(hsa[a[k]]), suma += calc(hsa[a[k]] ^= val[j]);
sumb -= calc(hsb[b[k]]), sumb += calc(hsb[b[k]] ^= val[j]);
}
mxl[i] = j;
}
while (q--) {
int l, r;
scanf("%d%d", &l, &r);
puts(l >= mxl[r] ? "Yes" : "No");
}
return 0;
}
CF799F Beautiful fountains rows
给定一个 \(n \times m\) 的矩阵 \(a\) ,其中第 \(i\) 行的 \(l_i \sim r_i\) 列为 \(1\) ,其余为 \(0\) 。
定义 \((a, b)\) 是好的当且仅当:
- \(\sum_{i = 1}^n \sum_{j = a}^b a_{i, j} > 0\) 。
- 对于任意第 \(i\) 行,\(\sum_{j = a}^b a_{i, j}\) 为 \(0\) 或奇数。
求所有好点对 \((a, b)\) 的 \(b - a + 1\) 之和。
\(n, m \leq 2 \times 10^5\)
看到奇偶性想到哈希。考虑给每一行的 \(1\) 映射到随机权值,那么一行和为奇数当且仅当区间异或和非 \(0\) 。
但是这样发现和判断 \(0\) 的条件不符,而且不好做到每一行都限制。考虑将奇数的限制也转化为区间异或和为 \(0\) ,即少统计一个有值位置。
由于行内 \(1\) 连续,考虑对于第 \(i\) 行令 \(l_i\) 处的哈希值为 \(0\) ,则两个条件的判定均变为 \(a + 1 \sim b\) 的哈希异或和为 \(0\) 。
那么对于整个矩阵,判定条件可以变为每一行的哈希异或和的异或和为 \(0\) ,再去掉全 \(0\) 的情况即可。
用 map
记录每个前缀的异或和可以做到 \(O(m \log m)\) 。
#include <bits/stdc++.h>
typedef unsigned long long ull;
typedef long long ll;
using namespace std;
const int N = 2e5 + 7;
ull val[N];
int c[N];
mt19937_64 myrand(time(0));
int n, m;
signed main() {
scanf("%d%d", &n, &m);
for (int i = 1; i <= n; ++i) {
int l, r;
scanf("%d%d", &l, &r);
ull w = myrand();
val[l + 1] ^= w, val[r + 1] ^= w;
++c[l], --c[r + 1];
}
for (int i = 1; i <= m; ++i)
val[i] ^= val[i - 1], c[i] += c[i - 1];
for (int i = 1; i <= m; ++i)
val[i] ^= val[i - 1];
map<ull, pair<int, ll> > mp;
ll ans = 0;
for (int i = 1; i <= m; ++i) {
auto &it = mp[val[i]];
++it.first, it.second += i - 1;
ans += 1ll * it.first * i - it.second;
}
for (int i = 1, len = 0; i <= m; ++i) {
len = (c[i] ? 0 : len + 1);
ans -= 1ll * len * (len + 1) / 2;
}
printf("%lld", ans);
return 0;
}
P6688 可重集
给定序列 \(a_{1 \sim n}\) ,\(m\) 次操作:
0 x y
:修改 \(a_x\) 为 \(y\) 。1 l r x y
:判断 \(a_{l \sim r}\) 与 \(a_{x \sim y}\) 是否本质相同,保证 \(r - l = y - x\) 。本质相同定义为将两序列升序排序后对应位置元素差相等。\(n, m \leq 10^6\)
若两序列本质相同,可以得到差值即为最小值的差。
设计可重集的哈希函数为 \(Hash(S) = \sum_{x \in S} base^x\) ,优势是无序、整体加减比较容易。
直接线段树维护即可做到 \(O(q \log n)\) 。
#include <bits/stdc++.h>
using namespace std;
const int Mod = 1e9 + 7;
const int N = 1e6 + 7;
int a[N];
mt19937 myrand(time(0));
int n, m, base;
inline int add(int x, int y) {
x += y;
if (x >= Mod)
x -= Mod;
return x;
}
inline int dec(int x, int y) {
x -= y;
if (x < 0)
x += Mod;
return x;
}
inline int mi(int a, int b) {
int res = 1;
for (; b; b >>= 1, a = 1ll * a * a % Mod)
if (b & 1)
res = 1ll * res * a % Mod;
return res;
}
namespace SMT {
pair<int, int> s[N << 2];
inline int ls(int x) {
return x << 1;
}
inline int rs(int x) {
return x << 1 | 1;
}
inline pair<int, int> operator + (const pair<int, int> &a, const pair<int, int> &b) {
return make_pair(min(a.first, b.first), add(a.second, b.second));
}
inline void pushup(int x) {
s[x] = s[ls(x)] + s[rs(x)];
}
void build(int x, int l, int r) {
if (l == r) {
s[x] = make_pair(a[l], mi(base, a[l]));
return;
}
int mid = (l + r) >> 1;
build(ls(x), l, mid), build(rs(x), mid + 1, r);
pushup(x);
}
void update(int x, int nl, int nr, int pos, int k) {
if (nl == nr) {
s[x] = make_pair(k, mi(base, k));
return;
}
int mid = (nl + nr) >> 1;
if (pos <= mid)
update(ls(x), nl, mid, pos, k);
else
update(rs(x), mid + 1, nr, pos, k);
pushup(x);
}
pair<int, int> query(int x, int nl, int nr, int l, int r) {
if (l <= nl && nr <= r)
return s[x];
int mid = (nl + nr) >> 1;
if (r <= mid)
return query(ls(x), nl, mid, l, r);
else if (l > mid)
return query(rs(x), mid + 1, nr, l, r);
else
return query(ls(x), nl, mid, l, r) + query(rs(x), mid + 1, nr, l, r);
}
} // namespace SMT
signed main() {
scanf("%d%d", &n, &m);
base = myrand() % (Mod - 1) + 1;
for (int i = 1; i <= n; ++i)
scanf("%d", a + i);
SMT::build(1, 1, n);
while (m--) {
int op;
scanf("%d", &op);
if (op) {
int l1, r1, l2, r2;
scanf("%d%d%d%d", &l1, &r1, &l2, &r2);
auto res1 = SMT::query(1, 1, n, l1, r1), res2 = SMT::query(1, 1, n, l2, r2);
if (res1.first > res2.first)
swap(res1, res2);
puts(1ll * res1.second * mi(base, res2.first - res1.first) % Mod == res2.second ? "YES" : "NO");
} else {
int x, k;
scanf("%d%d", &x, &k);
SMT::update(1, 1, n, x, k);
}
}
return 0;
}
P5278 算术天才⑨与等差数列
给定序列 \(a_{1 \sim n}\) ,\(m\) 次操作:
1 x y
:修改 \(a_x\) 为 \(y\) 。2 l r k
:判断 \(a_{l \sim r}\) 能否重排为公差为 \(k\) 的等差数列。\(n, m \leq 3 \times 10^5\)
可以重排时一般会将其视作集合,考虑套用可重集的哈希函数判断。
构造函数 calc(n, d)
表示长度为 \(n\) 、首项为 \(0\) 、公差为 \(d\) 的等差数列的哈希值,可以类似快速幂递归计算。
注意首项的平移,时间复杂度 \(O(q \log n)\) 。
#include <bits/stdc++.h>
using namespace std;
const int base = 233, Mod = 1e9 + 7;
const int N = 3e5 + 7;
int a[N];
int n, m;
inline int add(int x, int y) {
x += y;
if (x >= Mod)
x -= Mod;
return x;
}
inline int dec(int x, int y) {
x -= y;
if (x < 0)
x += Mod;
return x;
}
inline int mi(int a, int b) {
int res = 1;
for (; b; b >>= 1, a = 1ll * a * a % Mod)
if (b & 1)
res = 1ll * res * a % Mod;
return res;
}
int calc(int n, int d) { // a[1 ~ n], a[1] = 0
if (n == 1)
return 1;
else if (n & 1)
return add(1ll * calc(n - 1, d) * mi(base, d) % Mod, 1);
else
return 1ll * calc(n / 2, d) * add(mi(base, 1ll * d * (n / 2) % (Mod - 1)), 1) % Mod;
}
namespace SMT {
pair<int, int> s[N << 2];
inline int ls(int x) {
return x << 1;
}
inline int rs(int x) {
return x << 1 | 1;
}
inline pair<int, int> operator + (const pair<int, int> &a, const pair<int, int> &b) {
return make_pair(min(a.first, b.first), add(a.second, b.second));
}
inline void pushup(int x) {
s[x] = s[ls(x)] + s[rs(x)];
}
void build(int x, int l, int r) {
if (l == r) {
s[x] = make_pair(a[l], mi(base, a[l]));
return;
}
int mid = (l + r) >> 1;
build(ls(x), l, mid), build(rs(x), mid + 1, r);
pushup(x);
}
void update(int x, int nl, int nr, int pos, int k) {
if (nl == nr) {
s[x] = make_pair(k, mi(base, k));
return;
}
int mid = (nl + nr) >> 1;
if (pos <= mid)
update(ls(x), nl, mid, pos, k);
else
update(rs(x), mid + 1, nr, pos, k);
pushup(x);
}
pair<int, int> query(int x, int nl, int nr, int l, int r) {
if (l <= nl && nr <= r)
return s[x];
int mid = (nl + nr) >> 1;
if (r <= mid)
return query(ls(x), nl, mid, l, r);
else if (l > mid)
return query(rs(x), mid + 1, nr, l, r);
else
return query(ls(x), nl, mid, l, r) + query(rs(x), mid + 1, nr, l, r);
}
} // namespace SMT
signed main() {
scanf("%d%d", &n, &m);
for (int i = 1; i <= n; ++i)
scanf("%d", a + i);
SMT::build(1, 1, n);
int cnt = 0;
while (m--) {
int op;
scanf("%d", &op);
if (op == 1) {
int x, k;
scanf("%d%d", &x, &k);
x ^= cnt, k ^= cnt;
SMT::update(1, 1, n, x, k);
} else {
int l, r, k;
scanf("%d%d%d", &l, &r, &k);
l ^= cnt, r ^= cnt, k ^= cnt;
auto res = SMT::query(1, 1, n, l, r);
if (res.second == 1ll * calc(r - l + 1, k) * mi(base, res.first) % Mod)
++cnt, puts("Yes");
else
puts("No");
}
}
return 0;
}
CF1476E Kazaee
给定 \(a_{1 \sim n}\) ,\(m\) 次操作:
- 修改 \(a_x\) 为 \(k\) 。
- 查询 \(a_{l \sim r}\) 中所有数的出现次数为 \(k\) 的倍数。
\(n, m \leq 3 \times 10^5\)
这个查询一看就很难处理,但是可以发现的是这个查询和数字的具体值没有关系。
考虑随机将一个数随机映射为另一个数,注意到区间和是 \(k\) 的倍数是出现次数均为 \(k\) 的倍数的必要条件。
如果直接做可以随便卡,但是随机映射后就可以获得很好的正确性。
#include <bits/stdc++.h>
typedef long long ll;
using namespace std;
const int N = 3e5 + 7, K = 31;
map<int, int> mp;
ll id[N << 1][K];
int a[N];
mt19937 myrand(time(0));
int n, q, tot;
struct BIT {
ll c[N];
inline void update(int x, ll k) {
for (; x <= n; x += x & -x)
c[x] += k;
}
inline ll query(int x) {
ll res = 0;
for (; x; x -= x & -x)
res += c[x];
return res;
}
} bit[K];
signed main() {
scanf("%d%d", &n, &q);
for (int i = 1; i <= n; ++i) {
scanf("%d", a + i);
if (mp.find(a[i]) == mp.end()) {
mp[a[i]] = ++tot;
for (int j = 0; j < K; ++j)
id[tot][j] = myrand();
}
a[i] = mp[a[i]];
for (int j = 0; j < K; ++j)
bit[j].update(i, id[a[i]][j]);
}
while (q--) {
int op;
scanf("%d", &op);
if (op == 1) {
int x, k;
scanf("%d%d", &x, &k);
if (mp.find(k) == mp.end()) {
mp[k] = ++tot;
for (int i = 0; i < K; ++i)
id[tot][i] = myrand();
}
k = mp[k];
for (int i = 0; i < K; ++i)
bit[i].update(x, id[k][i] - id[a[x]][i]);
a[x] = k;
} else {
int l, r, k;
scanf("%d%d%d", &l, &r, &k);
bool flag = true;
for (int i = 0; i < K; ++i)
if ((bit[i].query(r) - bit[i].query(l - 1)) % k) {
flag = false;
break;
}
puts(flag ? "YES" : "NO");
}
}
return 0;
}
QOJ6504. Flower's Land 2
一个序列是合法的,当且仅当可以重复如下操作删空:删去两个相邻相同的元素。
给出长度为 \(n\) 的序列,值域为 \(\{ 0, 1, 2 \}\) 。\(q\) 次操作,操作有:
- 给出 \(l, r\) ,对于所有 \(l \leq i \leq r\) ,令 \(a_i \to (a_i + 1) \bmod 3\) 。
- 给出 \(l, r\) ,判定 \(a_{l \sim r}\) 是否合法。
\(n, q \leq 5 \times 10^5\)
先考虑如何判定一个区间合法,显然可以用栈维护消除的过程,但是这个过程很难扩展以支持区间修改。
考虑将 \(0, 1, 2\) 随机映射为 \(A_{0, 1, 2}\) ,由于每次消去的位置的奇偶性一定不同,因此令奇数位为 \(A_{0, 1, 2}\) ,偶数位为 \(A^{-1}_{0, 1, 2}\) ,若区间连乘为 \(1\) 则认为该区间合法。
但是映射到数字会产生一个问题,因为数字的连乘是具有交换律的,因此其判定不满足相邻的限制。
自然想到映射到随机矩阵即可,线段树维护当前区间 \(+0, +1, +2\) 后的区间乘积即可,时间复杂度 \(O((n + q \log n) B^3)\) ,其中 \(B\) 为矩阵边长。
#include <bits/stdc++.h>
using namespace std;
const int Mod = 1e9 + 7;
const int N = 5e5 + 7, B = 2;
char str[N];
mt19937 myrand(time(0));
int n, q;
inline int add(int x, int y) {
x += y;
if (x >= Mod)
x -= Mod;
return x;
}
inline int dec(int x, int y) {
x -= y;
if (x < 0)
x += Mod;
return x;
}
inline int mi(int a, int b) {
int res = 1;
for (; b; b >>= 1, a = 1ll * a * a % Mod)
if (b & 1)
res = 1ll * res * a % Mod;
return res;
}
struct Matrix {
int a[B][B];
inline Matrix() {
memset(a, 0, sizeof(a));
}
inline Matrix operator * (const Matrix &rhs) {
Matrix res;
for (int i = 0; i < B; ++i)
for (int j = 0; j < B; ++j)
for (int k = 0; k < B; ++k)
res.a[i][k] = add(res.a[i][k], 1ll * a[i][j] * rhs.a[j][k] % Mod);
return res;
}
} a[2][3];
inline Matrix Gauss(Matrix mt) {
static int g[B][B * 2];
for (int i = 0; i < B; ++i) {
memcpy(g[i], mt.a[i], sizeof(int) * B);
memset(g[i] + B, 0, sizeof(int) * B);
}
for (int i = 0; i < B; ++i)
g[i][i + B] = 1;
for (int i = 0; i < B; ++i) {
int inv = mi(g[i][i], Mod - 2);
for (int j = 0; j < B; ++j)
if (j != i) {
int div = 1ll * g[j][i] * inv % Mod;
for (int k = i; k < B * 2; ++k)
g[j][k] = dec(g[j][k], 1ll * div * g[i][k] % Mod);
}
}
for (int i = 0; i < B; ++i) {
int inv = mi(g[i][i], Mod - 2);
for (int j = 0; j < B; ++j)
mt.a[i][j] = 1ll * g[i][B + j] * inv % Mod;
}
return mt;
}
namespace SMT {
Matrix mt[N << 2][3];
int tag[N << 2];
inline int ls(int x) {
return x << 1;
}
inline int rs(int x) {
return x << 1 | 1;
}
inline void pushup(int x) {
for (int i = 0; i < 3; ++i)
mt[x][i] = mt[ls(x)][i] * mt[rs(x)][i];
}
inline void spread(int x, int k) {
rotate(mt[x], mt[x] + (k % 3), mt[x] + 3), tag[x] = (tag[x] + k) % 3;
}
inline void pushdown(int x) {
if (tag[x])
spread(ls(x), tag[x]), spread(rs(x), tag[x]), tag[x] = 0;
}
void build(int x, int l, int r) {
if (l == r) {
for (int i = 0; i < 3; ++i)
mt[x][i] = a[l & 1][((str[l] & 15) + i) % 3];
return;
}
int mid = (l + r) >> 1;
build(ls(x), l, mid), build(rs(x), mid + 1, r);
pushup(x);
}
void update(int x, int nl, int nr, int l, int r, int k) {
if (l <= nl && nr <= r) {
spread(x, k);
return;
}
pushdown(x);
int mid = (nl + nr) >> 1;
if (l <= mid)
update(ls(x), nl, mid, l, r, k);
if (r > mid)
update(rs(x), mid + 1, nr, l, r, k);
pushup(x);
}
Matrix query(int x, int nl, int nr, int l, int r) {
if (l <= nl && nr <= r)
return mt[x][0];
pushdown(x);
int mid = (nl + nr) >> 1;
if (r <= mid)
return query(ls(x), nl, mid, l, r);
else if (l > mid)
return query(rs(x), mid + 1, nr, l, r);
else
return query(ls(x), nl, mid, l, r) * query(rs(x), mid + 1, nr, l, r);
}
} // namespace SMT
signed main() {
for (int i = 0; i < 3; ++i) {
for (int j = 0; j < B; ++j)
for (int k = 0; k < B; ++k)
a[0][i].a[j][k] = myrand() % Mod;
a[1][i] = Gauss(a[0][i]);
}
scanf("%d%d%s", &n, &q, str + 1);
SMT::build(1, 1, n);
while (q--) {
int op, l, r;
scanf("%d%d%d", &op, &l, &r);
if (op == 1)
SMT::update(1, 1, n, l, r, 1);
else {
Matrix ans = SMT::query(1, 1, n, l, r);
bool flag = true;
for (int i = 0; i < B && flag; ++i)
for (int j = 0; j < B && flag; ++j)
flag &= (ans.a[i][j] == (i == j));
puts(flag ? "Yes" : "No");
}
}
return 0;
}
CF1641D Two Arrays
给定 \(n\) 个长度为 \(m\) 的序列,每个序列有各自的价值。找到一对元素不交的序列,最小化价值和。
\(n \leq 10^5\) ,\(m \leq 5\)
考虑将所有数字随机映射到 \([0, B)\) ,然后直接做高维前缀和可以做到 \(O(nm + B \times 2^B)\) 。
取 \(B = 20\) ,随机 \(100\) 次正确率就很高了。
#include <bits/stdc++.h>
using namespace std;
const int inf = 0x3f3f3f3f;
const int N = 1e5 + 7, M = 5, B = 18;
vector<int> vec;
int a[N][M], val[N], f[1 << B];
mt19937 myrand(time(0));
int n, m;
signed main() {
scanf("%d%d", &n, &m);
for (int i = 1; i <= n; ++i) {
for (int j = 0; j < m; ++j)
scanf("%d", a[i] + j), vec.emplace_back(a[i][j]);
scanf("%d", val + i);
}
sort(vec.begin(), vec.end());
vec.erase(unique(vec.begin(), vec.end()), vec.end());
for (int i = 1; i <= n; ++i)
for (int j = 0; j < m; ++j)
a[i][j] = lower_bound(vec.begin(), vec.end(), a[i][j]) - vec.begin();
vector<int> id(vec.size());
iota(id.begin(), id.end(), 0);
int ans = 2e9 + 1;
for (int times = 1; times <= 150; ++times) {
shuffle(id.begin(), id.end(), myrand);
memset(f, inf, sizeof(f));
for (int i = 1; i <= n; ++i) {
int s = 0;
for (int j = 0; j < m; ++j)
s |= 1 << (id[a[i][j]] % B);
if (__builtin_popcount(s) == m)
f[s] = min(f[s], val[i]);
}
for (int i = 0; i < B; ++i)
for (int j = 0; j < (1 << B); ++j)
if (j >> i & 1)
f[j] = min(f[j], f[j ^ (1 << i)]);
for (int i = 0; i < (1 << B); ++i)
if (f[i] != inf && f[(1 << B) - 1 - i] != inf)
ans = min(ans, f[i] + f[(1 << B) - 1 - i]);
}
printf("%d", ans == 2e9 + 1 ? -1 : ans);
return 0;
}
CF2003F Turtle and Three Sequences
有 \(n\) 个三元组 \((a_i, b_i, c_i)\) ,选 \(m\) 个三元组满足 \(a\) 单调不降且 \(b\) 互不相同,最大化 \(\sum c_i\) 或报告无解。
\(n \leq 3000\) ,\(m \le 5\)
发现 \(b\) 互不相同不好处理,考虑将 \(b\) 随机映射到 \(0 \sim m - 1\) 上的整数,则可以状压 DP。
一次映射的成功率为 \(\frac{m!}{m^m}\) ,期望随机 \(\frac{m^m}{m!}\) 次即可,单次 DP 可以树状数组优化做到 \(O(2^m n \log n)\) 。
#include <bits/stdc++.h>
using namespace std;
const int inf = 0x3f3f3f3f;
const int N = 3e3 + 7, M = 5;
int a[N], b[N], c[N], hs[N];
mt19937 myrand(time(0));
int n, m;
struct BIT {
int c[N];
inline void update(int x, int k) {
for (; x <= n; x += x & -x)
c[x] = max(c[x], k);
}
inline int query(int x) {
int res = -inf;
for (; x; x -= x & -x)
res = max(res, c[x]);
return res;
}
} bit[1 << M];
signed main() {
scanf("%d%d", &n, &m);
for (int i = 1; i <= n; ++i)
scanf("%d", a + i);
for (int i = 1; i <= n; ++i)
scanf("%d", b + i);
for (int i = 1; i <= n; ++i)
scanf("%d", c + i);
int ans = -1;
while (clock() / CLOCKS_PER_SEC < 2.8) {
for (int i = 1; i <= n; ++i)
hs[i] = myrand() % m;
for (int s = 0; s < (1 << m); ++s)
memset(bit[s].c, -0x3f, sizeof(bit[s].c));
bit[0].update(1, 0);
for (int i = 1; i <= n; ++i)
for (int s = 0; s < (1 << m); ++s)
if (~s >> hs[b[i]] & 1)
bit[s | (1 << hs[b[i]])].update(a[i], bit[s].query(a[i]) + c[i]);
ans = max(ans, bit[(1 << m) - 1].query(n));
}
printf("%d", ans);
return 0;
}
利用随机化性质
- 随机父节点树的期望树高为 \(O(\log n)\) 。
- 随机树的期望直径为 \(O(\sqrt{n})\) 。
- 对于随机排列 \(p_{1 \sim n}\) ,\(p_i\) 作为前缀最大值的概率为 \(\frac{1}{i}\) ,前缀最大值的期望个数为 \(O(\sum \frac{1}{i}) = O(\ln n)\) 。
- 在 \([1, n]\) 中随机选数,选出两个相同数的期望轮数近似 \(2 \sqrt{n}\) 。
- 对于随机序列 \(x_{1 \sim n}\) ,其中 \(x_i\) 等概率从 \(\{ -1, 1 \}\) 中选取,则前缀和的最大值为 \(O(\sqrt{n})\) 级别。
P9047 [PA 2021] Poborcy podatkowi
给定一棵 \(n\) 个点的树,可以选择若干条长度为 \(4\) 的边不相交链(可以不选),最大化选出边集的边权和。
\(n \leq 2 \times 10^5\)
设 \(f_{u, i}\) 表示 \(u\) 子树内留一条 \(u\) 开端的长度为 \(i\) 的链的答案,转移时对儿子 \(f_{v, i}\) 分类讨论:
- \(i = 3\) :直接加上 \((u, v)\) 尝试贡献答案即可。
- \(i = 0 / 2\) :需要 \(i = 0\) 与 \(i = 2\) 两两匹配。
- \(i = 1\) :内部随意两两匹配。
因此考虑设 \(g_{i, j, k}\) 表示考虑了前 \(i\) 个儿子,\(j\) 记录未匹配的 \(1\) 的数量,\(k\) (带正负,记 \(0\) 为 \(1\) ,\(2\) 为 \(-1\) )记录未匹配的 \(0 / 2\) 的数量。由于 \(j \in \{ 0, 1 \}\) ,因此不难暴力转移做到 \(O(n^2)\) 。
考虑优化,对于一组儿子 \((0, 2)\) 匹配或不匹配的方案,可以将其写成一个 \(\in \{-1, 0, 1 \}\) 的序列,将其随机打乱后,最大前缀和是 \(O(\sqrt{|son(u)|})\) 级别的。因此考虑随机打乱儿子,这样最优方案的最后一维转移上下界仅有 \(O(\sqrt{|son(u)|})\) 级别,因而可以忽略范围外的 DP 值。
将最后一维的阈值设为 \(O(\sqrt{|son(u)|})\) 级别即可做到时间复杂度 \(O(n \sqrt{n})\) 。
#include <bits/stdc++.h>
typedef long long ll;
using namespace std;
const ll inf = 1e18;
const int N = 2e5 + 7;
struct Graph {
vector<pair<int, int> > e[N];
inline void insert(int u, int v, int w) {
e[u].emplace_back(v, w);
}
} G;
ll f[N][4];
mt19937 myrand(time(0));
int n;
void dfs(int u, int fa) {
for (auto it : G.e[u]) {
int v = it.first;
if (v != fa)
dfs(v, u);
}
shuffle(G.e[u].begin(), G.e[u].end(), myrand);
int B = sqrt(G.e[u].size()) * 2;
vector<ll> g[2] = {vector<ll>(B << 1 | 1, -inf), vector<ll>(B << 1 | 1, -inf)};
g[0][0 + B] = 0;
for (auto it : G.e[u]) {
int v = it.first, w = it.second;
if (v == fa)
continue;
vector<ll> h[2] = {vector<ll>(B << 1 | 1, -inf), vector<ll>(B << 1 | 1, -inf)};
for (int i = 0; i <= 1; ++i)
for (int j = -B; j < B; ++j) {
h[i][j + B] = max(h[i][j + B], g[i][j + B] + f[v][0]);
for (int k = 0; k <= 3; ++k) {
int d = (k ? (k == 2 ? -1 : 0) : 1);
if (-B <= j + d && j + d <= B)
h[i ^ (k == 1)][j + d + B] = max(h[i ^ (k == 1)][j + d + B], g[i][j + B] + f[v][k] + w);
}
}
swap(g, h);
}
f[u][0] = g[0][0 + B], f[u][1] = g[0][1 + B], f[u][2] = g[1][0 + B], f[u][3] = g[0][-1 + B];
}
signed main() {
scanf("%d", &n);
for (int i = 1; i < n; ++i) {
int u, v, w;
scanf("%d%d%d", &u, &v, &w);
G.insert(u, v, w), G.insert(v, u, w);
}
dfs(1, 0);
printf("%lld\n", f[1][0]);
return 0;
}
异或
二维平面上有 \(n\) 个关键点,每个点有一个非负整数权值 \(a_i\) 。
定义一个点 \(P\) 的最小半径 \(f(P)\) 为满足以下条件的 \(d\) 的最小值:
- 存在一个关键点集合,满足每个关键点到 \(P\) 的距离 \(\leq d\) ,且该关键点集合内权值的异或和 \(\geq k\) ,其中 \(k\) 为给定的常数。
- 若不存在合法的 \(d\) ,则 \(f(P) = +\infty\) 。
求二维平面上所有点(不一定是关键点)的最小半径的最小值。
\(n \leq 3000\) ,\(0 \leq a_i, k < 2^{20}\)
显然最优解是一个包含若干点的最小点覆盖。考虑枚举一个点 \(P\) 在关键点集合内,二分半径 \(r\) ,判定是否存在一个包含 \(P\) 、半径为 \(r\) 的合法的圆。
不妨钦定 \(P\) 在圆周上,这样可能的圆心构成一个以 \(P\) 为圆心、半径为 \(r\) 的圆。对于点 \(Q\) ,若 \(|PQ| > r\) ,则 \(Q\) 不可能被统计到,否则 \(Q\) 对应的合法的圆心范围是一段圆弧。
判定可以用线性基维护,这里无需用线段树分治,只要用记录每个数的删除时间,用前缀线性基维护即可。
直接做是 \(O(n \log V \times n(\log n + \log A))\) 的,无法通过。
考虑经典结论:随机排列的前缀最大值只有 \(O(\log n)\) 个。考虑将 \(n\) 个点随机 shuffle
,记当前答案为 \(ans\) ,遍历到第 \(i\) 个点是先判定 \(r = ans\) 是否合法,若不合法则退出,否则继续二分。
不难发现期望二分次数是 \(O(\log n)\) 次的,时间复杂度 \(O((n + \log n \log V) \times n (\log n + \log A))\) ,可以通过。
#include <bits/stdc++.h>
using namespace std;
const double eps = 1e-9;
const int N = 3e3 + 7, B = 21;
struct Point {
double x, y;
inline friend double dist(const Point &a, const Point &b) {
return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}
} p[N];
struct LinearBasis {
int p[B], d[B];
inline LinearBasis() {
memset(p, 0, sizeof(p)), memset(d, -1, sizeof(d));
}
inline void insert(int x, int k) {
for (int i = B - 1; ~i; --i)
if (x >> i & 1) {
if (!p[i]) {
p[i] = x, d[i] = k;
return;
} else {
if (k > d[i])
swap(x, p[i]), swap(k, d[i]);
x ^= p[i];
}
}
}
inline void remove(int k) {
for (int i = B - 1; ~i; --i)
if (p[i] && d[i] == k)
p[i] = 0, d[i] = -1;
}
inline int querymax() {
int res = 0;
for (int i = B - 1; ~i; --i)
if (p[i] && (~res >> i & 1))
res ^= p[i];
return res;
}
} lb;
int a[N];
mt19937 myrand(time(0));
int n, m;
inline bool check(int x, double r) {
vector<pair<double, int> > upd;
for (int i = 1; i <= n; ++i) {
if (p[i].x == p[x].x && p[i].y == p[x].y) {
upd.emplace_back(-1e9, i);
continue;
} else if (dist(p[i], p[x]) > r * 2)
continue;
double a = -2 * (p[i].x - p[x].x), b = -2 * (p[i].y - p[x].y),
c = p[i].x * p[i].x + p[i].y * p[i].y - p[x].x * p[x].x - p[x].y * p[x].y,
K = -a / b, B = -c / b;
auto solve = [](double a, double b, double c) {
double delta = b * b - 4 * a * c;
return make_pair((-b - sqrt(delta)) / 2 / a, (-b + sqrt(delta)) / 2 / a);
};
pair<double, double> X = solve(1 + K * K, -2 * p[x].x + 2 * (B - p[x].y) * K,
p[x].x * p[x].x + (B - p[x].y) * (B - p[x].y) - r * r),
Y = make_pair(K * X.first + B, K * X.second + B);
auto angle = [&](double px, double py) {
return atan2(py - p[x].y, px - p[x].x);
};
double l = angle(X.first, Y.first), mid = angle(p[i].x, p[i].y), r = angle(X.second, Y.second);
if (l > r)
swap(l, r);
if (l <= mid && mid <= r) {
upd.emplace_back(l, i), upd.emplace_back(r, -i);
upd.emplace_back(l + 1e9, i), upd.emplace_back(r + 1e9, -i);
} else
upd.emplace_back(r, i), upd.emplace_back(l + 1e9, -i), upd.emplace_back(r + 1e9, i);
}
sort(upd.begin(), upd.end(), [](const pair<double, int> &a, const pair<double, int> &b) {
return a.first == b.first ? a.second > b.second : a.first < b.first;
});
vector<int> lst(n + 1, upd.size() + 1), tim(upd.size());
for (int i = upd.size() - 1; ~i; --i) {
int x = upd[i].second;
if (x > 0)
tim[i] = lst[x], lst[x] = upd.size() + 1;
else
lst[-x] = i;
}
lb = LinearBasis();
for (int i = 0; i < upd.size(); ++i) {
int x = upd[i].second;
if (x > 0)
lb.insert(a[x], tim[i]);
else
lb.remove(i);
if (lb.querymax() >= m)
return true;
}
return false;
}
signed main() {
scanf("%d%d", &n, &m);
for (int i = 1; i <= n; ++i)
scanf("%lf%lf%d", &p[i].x, &p[i].y, a + i), lb.insert(a[i], 0);
if (lb.querymax() < m)
return puts("-1"), 0;
vector<int> id(n);
iota(id.begin(), id.end(), 1), shuffle(id.begin(), id.end(), myrand);
double ans = 1e9;
for (int x : id) {
if (!check(x, ans))
continue;
double l = 0, r = ans;
while (r - l > eps) {
double mid = (l + r) / 2;
if (check(x, mid))
ans = r = mid;
else
l = mid;
}
}
printf("%.9lf", ans);
return 0;
}
P8108 [Cnoi2021] 绀珠传说
给定一个 \(n \times n\) 的网格,其中 \(1 \sim n\) 的每个整数都恰好出现 \(n\) 次,保证矩阵随机生成。
定义一次操作为:选择最底行的一个连续段将其删去,并将其上方的部分下移。
求删完所有数的最少操作次数。
\(n \leq 1000\)
考虑将最优方案中一起被删去的点连成一个连通块,对于相邻的两列之间的连边,限制是不能交叉,因此答案即为 \(n^2\) 减去所有相邻列的 LCS。
由于数据随机,因此相邻列的相同元素对的数量为 \(O(n)\) 级别(每个数期望出现 \(O(1)\) 次),树状数组优化 DP 即可做到 \(O(n^2 \log n)\) 。
#include <bits/stdc++.h>
using namespace std;
const int N = 1e3 + 7;
int a[N][N];
int n;
namespace BIT {
int c[N];
inline void update(int x, int k) {
for (; x <= n; x += x & -x)
c[x] = max(c[x], k);
}
inline int query(int x) {
int res = 0;
for (; x; x -= x & -x)
res = max(res, c[x]);
return res;
}
} // namespace BIT
signed main() {
scanf("%d", &n);
for (int i = 1; i <= n; ++i)
for (int j = 1; j <= n; ++j)
scanf("%d", a[i] + j);
int ans = n * n;
for (int i = 1; i < n; ++i) {
vector<vector<int> > vec(n + 1);
for (int j = 1; j <= n; ++j)
vec[a[j][i]].emplace_back(j);
vector<pair<int, int> > nd;
for (int j = 1; j <= n; ++j)
for (int it : vec[a[j][i + 1]])
nd.emplace_back(it, j);
sort(nd.begin(), nd.end(), [](const pair<int, int> a, const pair<int, int> b) {
return a.first == b.first ? a.second > b.second : a.first < b.first;
});
memset(BIT::c + 1, 0, sizeof(int) * n);
for (auto it : nd)
BIT::update(it.second, BIT::query(it.second - 1) + 1);
ans -= BIT::query(n);
}
printf("%d", ans);
return 0;
}