ACM模板整理

文章目录
  1. 1. 姿势
    1. 1.1. 手动扩栈
    2. 1.2. Code::Blocks输入重定向
  2. 2. 计算几何
  3. 3. 读入输出优化
    1. 3.1. 简单版
    2. 3.2. 复杂版
  4. 4. KMP
  5. 5. Hungary
  6. 6. 逆元
    1. 6.1. 扩展欧几里得
    2. 6.2. 费马小定理
    3. 6.3. 线性递推
    4. 6.4. 线性求阶乘的逆元
  7. 7. 线性筛同时求质数、欧拉函数、莫比乌斯函数、莫比乌斯函数前缀和
  8. 8. \(1 \leq x \leq n\), \(1 \leq y \leq m\), \(gcd(x,y)=k\)的个数
  9. 9. FWT
  10. 10. pb_ds库
    1. 10.1. priority_queue
    2. 10.2. tree
    3. 10.3. hash_table
  11. 11. 最小费用最大流
  12. 12. k维最远曼哈顿距离
  13. 13. Java常用库
  14. 14. 对拍
    1. 14.1. check.bat (Windows BAT)
    2. 14.2. check.sh (Linux shell)
    3. 14.3. rand.cpp
  15. 15. FFT
  16. 16. 数值积分
  17. 17. Dijkstra+pb_ds
  18. 18. AC自动机
  19. 19. 球盒问题
    1. 19.1. 第二类斯特林数
    2. 19.2. 球、盒均相同的DP方法
  20. 20. 第一类斯特林数
  21. 21. Catalan数
  22. 22. 错排
  23. 23. 康拓展开
    1. 23.1. 逆运算
  24. 24. rope库
  25. 25. 莫比乌斯反演
    1. 25.1. 定理1
    2. 25.2. 定理2
  26. 26. 数位DP
  27. 27. 大整数类
  28. 28. 1~n中与m互质的数的和、平方和
  29. 29. 欧拉函数
  30. 30. \(ax+by=c\)的整数解
  31. 31. Lucas定理
  32. 32. 模线性方程组
    1. 32.1. \(m_i\)两两互质
    2. 32.2. \(m_i\)不满足两两互质
  33. 33. Splay
  34. 34. 除法取整分块
  35. 35. Heron法开方(其实也是Newton法)
  36. 36. Newton法开方
    1. 36.1. 整数
    2. 36.2. 小数
  37. 37. 平方剩余
    1. 37.1. 模数为奇质数
    2. 37.2. 模数为奇数
  38. 38. 线段树
  39. 39. 后缀自动机
  40. 40. Tarjan强连通分量
  41. 41. Tarjan割点与桥
  42. 42. 树分治
  43. 43. KD树
  44. 44. 树链剖分
  45. 45. 树状数组
  46. 46. 最大密度子图
  47. 47. 最大权闭合子图
  48. 48. 约瑟夫环
  49. 49. 博弈论
    1. 49.1. Bash游戏
    2. 49.2. Nim游戏
    3. 49.3. Wythoff游戏
  50. 50. 01分数规划
  51. 51. 表达式处理
    1. 51.1. 中缀表达式直接求值
    2. 51.2. 中缀表达式转后缀表达式
  52. 52. Miller-Rabin素性测试
  53. 53. std::regex

自用模板。

姿势

手动扩栈

g++

1
2
3
int size = 256 << 20; //256M
char* p = (char*)malloc(size) + size;
__asm__("movl %0, %%esp\n" :: "r"(p));

g++ x64

1
2
3
4
5
6
7
8
9
extern void main2() __asm__("main2");
void main2() {
...
exit(0);
}
int main() {
int size = 256 << 20;
__asm__ __volatile__("movq %0, %%rsp\njmp main2"::"r"((char*)malloc(size)+size));
}

c++

1
#pragma comment(linker, "/STACK:102400000,102400000")

Code::Blocks输入重定向

Windows下

Tools:

1
2
3
4
Name: redirect input
Executable: cmd.exe
Parameters: /C echo ========BUILDING======== && echo. && g++ -g -Wall -std=c++11 ${ACTIVE_EDITOR_FILENAME} -o ${TARGET_OUTPUT_BASENAME}.exe && echo ========START======== && echo. && ${TARGET_OUTPUT_BASENAME}.exe < in.txt
Working directory: ${PROJECT_DIR}

Linux下

Tools:

1
2
3
4
Name: redirect input
Executable: bash
Parameters: -c "echo ========BUILDING======== && g++ -g -Wall -std=c++11 ${ACTIVE_EDITOR_FILENAME} -o ${TARGET_OUTPUT_BASENAME}.o && echo ========START======== && ./${TARGET_OUTPUT_BASENAME}.o < in.txt"
Working directory: ${PROJECT_DIR}

设置快捷键: Settings -> Editor -> Keyboard shortcuts

计算几何

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <vector>
#include <cassert>
#include <utility>
namespace Geometry {
using namespace std;
typedef double real; // 需要时可以转为long double
const real EPS = 1e-6; // 精度控制
const real PI = acos((real)-1);
// 浮点数比较
// @param x: 任意实数
// @return 0: x==0, 1: x>0, -1: x<0
int dcmp(real x) {
return (x>EPS) - (x<-EPS);
}
int dcmp(real x, real y) {
return dcmp(x - y);
}
struct Point {
real x, y, z;
int id;
Point(real x = 0, real y = 0, real z = 0) :x(x), y(y), z(z) {}
// 强制转换为real,表示点到原点距离
operator real() const {
return sqrt(x*x + y * y + z * z);
}
// 平面向量与x轴正向夹角
// @return 方向角,(-pi, pi]
real direction() const {
return atan2(y, x);
}
};
typedef Point Vector;
Vector operator+ (Vector a, Vector b) {
return Vector(a.x + b.x, a.y + b.y, a.z + b.z);
}
Vector operator- (Vector a, Vector b) {
return Vector(a.x - b.x, a.y - b.y, a.z - b.z);
}
Vector operator- (Vector a) {
return Vector(-a.x, -a.y, -a.z);
}
Vector operator* (Vector a, real b) {
return Vector(a.x*b, a.y*b, a.z*b);
}
Vector operator* (real a, Vector b) {
return Vector(a*b.x, a*b.y, a*b.z);
}
// 向量点积
real operator* (Vector a, Vector b) {
return a.x*b.x + a.y*b.y + a.z*b.z;
}
Vector operator/ (Vector a, real b) {
return Vector(a.x / b, a.y / b, a.z / b);
}
// 向量叉积
Vector cross(Vector a, Vector b) {
return Vector(a.y*b.z - a.z*b.y, a.z*b.x - a.x*b.z, a.x*b.y - a.y*b.x);
}
// 平面向量叉积的值
real Cross(Vector a, Vector b) {
return a.x*b.y - a.y*b.x;
}
// 向量混合积=axb*c
real mix(Vector a, Vector b, Vector c) {
return cross(a, b)*c;
}
// 绕z轴旋转向量
Vector rotate(Vector v, real theta) {
return Vector(v.x*cos(theta) - v.y*sin(theta), v.x*sin(theta) + v.y*cos(theta), v.z);
}
real distance(Point a, Point b) {
return (real)(a - b);
}
// 平面向量v1->v2的旋转角
// @return (-pi, pi]
real rotAngle(Vector v1, Vector v2) {
real rst = v2.direction() - v1.direction();
if (rst>PI) rst -= PI * 2;
if (rst <= -PI) rst += PI * 2;
return rst;
}
// 向量夹角
// @return [0, pi]
real angle(Vector v1, Vector v2) {
real costheta = v1 * v2 / (real)v1 / (real)v2;
return acos(costheta);
}
// 方向和长度生成向量
Vector zoom(Vector v, real length) {
return v * length / (real)v;
}
// 有向直线,方向为a->b
struct Line {
Point a, b;
Line(Point a = Point(), Point b = Point()) :a(a), b(b) {}
// 方向向量
Vector direction() const {
return b - a;
}
};
// 点到直线距离
real distance(Point p, Line l) {
return abs((real)cross(l.b - l.a, p - l.a)) / real(l.b - l.a);
}
// 两直线位置关系
// @return 0:平行 1:相交 2:异面
int position(Line l1, Line l2) {
Vector v1 = l1.direction(), v2 = l2.direction();
if (dcmp((real)cross(v1, v2)) == 0) return 0;
if (dcmp(mix(l2.a - l1.a, v1, v2)) == 0) return 1;
return 2;
}
// 两直线交点
Point intersect(Line l1, Line l2) {
assert(position(l1, l2) == 1);
Vector c1 = cross(l1.direction(), l2.a - l1.a), c2 = cross(l1.direction(), l2.direction());
real sgn = dcmp(c1*c2)>0 ? -1 : 1;
return l2.a + l2.direction() * (real)c1 / (real)c2 * sgn;
}
// 判断三点共线
bool colinear(Point a, Point b, Point c) {
return dcmp((real)cross(b - a, c - a)) == 0;
}
// 平面点与直线的位置关系
// 0: 在直线上, 1: 点在直线左侧, -1: 点在直线右侧
int position(Point p, Line l) {
return dcmp(Cross(l.b - l.a, p - l.a));
}
// 常用坐标比较函数
bool cmpxyz(Point a, Point b) {
if (dcmp(a.x - b.x) != 0) return a.x<b.x;
else if (dcmp(a.y - b.y) != 0) return a.y<b.y;
else return dcmp(a.z - b.z)<0;
}

// 有向线段
struct Segment {
Point a, b;
Segment(Point a = Point(), Point b = Point()) : a(a), b(b) {}
// 方向向量
Vector direction() const {
return b - a;
}
// 线段长度
operator real() const {
return (real)direction();
}
};
// 点到线段距离
real distance(Point p, Segment seg) {
Vector ap = p - seg.a, bp = p - seg.b;
if (ap*(seg.b - seg.a) <= 0) return (real)ap;
if (bp*(seg.a - seg.b) <= 0) return (real)bp;
return distance(p, Line(seg.a, seg.b));
}
// 点和线段的位置关系
// @return true:在线段上 false:不在
bool position(Point p, Segment seg) {
return dcmp(distance(p, seg)) == 0;
}
// 线段相交,端点处相交也算
bool intersected(Segment a, Segment b) {
if (position(Line(a.a, a.b), Line(b.a, b.b)) != 1) return false;
return dcmp(cross(a.direction(), b.a - a.a)*cross(a.direction(), b.b - a.a)) <= 0 &&
dcmp(cross(b.direction(), a.a - b.a)*cross(b.direction(), a.b - b.a)) <= 0;
}

// 平面简单多边形
struct Polygon2D {
vector<Point> vtx;
// 按逆时针顺序给出顶点
Polygon2D(vector<Point> vertex = vector<Point>()) :vtx(vertex) {}
// 第i条边
// @param i: 0~n-1
Segment side(int i) const {
if (i == vtx.size() - 1) return Segment(vtx[vtx.size() - 1], vtx[0]);
return Segment(vtx[i], vtx[i + 1]);
}
// 面积
real area() const {
real rst = 0;
int sz = vtx.size();
for (int i = 0; i<sz; i++) {
rst += Cross(vtx[i], vtx[(i + 1) % sz]);
}
return rst / 2;
}
// 重心
Point cofg() const {
Point rst;
real ar = 0;
int sz = vtx.size();
for (int i = 0; i<sz; i++) {
real temp = Cross(vtx[i], vtx[(i + 1) % sz]);
rst = rst + (vtx[i] + vtx[(i + 1) % sz])*temp;
ar += temp;
}
return rst / ar / 3.;
}
// 周长
real circumference() const {
real rst = 0;
int sz = vtx.size();
for (int i = 0; i<sz; i++) {
rst += (real)side(i);
}
return rst;
}
// 凸包算法将点按逆时针排序
void arrange() {
sort(vtx.begin(), vtx.end(), cmpxyz);
vector<Point> p;
for (const Point& it : vtx) {
while (p.size() >= 2) {
Line last(p[p.size() - 2], p[p.size() - 1]);
if (position(it, last) == -1) {
p.pop_back();
}
else break;
}
p.push_back(it);
}
for (vector<Point>::const_reverse_iterator it = ++vtx.rbegin(); it != vtx.rend(); ++it) {
while (p.size() >= 2) {
Line last(p[p.size() - 2], p[p.size() - 1]);
if (position(*it, last) == -1) {
p.pop_back();
}
else break;
}
p.push_back(*it);
}
p.pop_back();
vtx = move(p);
}
};
// 点与多边形的位置关系
// @return -1:内 0:上 1:外
int position(Point p, Polygon2D c) {
int n = c.vtx.size();
int cnt = 0;
for (int i = 0; i < n; i++) {
Segment seg = c.side(i);
if (position(p, seg)) return 0;
int k = dcmp(Cross(seg.direction(), p - seg.a));
int d1 = dcmp(seg.a.y, p.y);
int d2 = dcmp(seg.b.y, p.y);
if (k>0 && d1 <= 0 && d2>0) cnt++;
if (k<0 && d2 <= 0 && d1>0) cnt--;
}
if (cnt) return -1;
else return 1;
}

// 平面圆
struct Circle2D {
Point ct;
real r;
Circle2D(Point center = Point(), real radius = 0) :ct(center), r(radius) {}
// 过三点的圆
Circle2D(Point a, Point b, Point c) {
real x1 = a.x, y1 = a.y, x2 = b.x, y2 = b.y, x3 = c.x, y3 = c.y;
real a11 = 2 * (x3 - x2);
real a12 = 2 * (y3 - y2);
real a21 = 2 * (x2 - x1);
real a22 = 2 * (y2 - y1);
real b1 = x3 * x3 - x2 * x2 + y3 * y3 - y2 * y2;
real b2 = x2 * x2 - x1 * x1 + y2 * y2 - y1 * y1;
real d = a11 * a22 - a12 * a21;
real d1 = b1 * a22 - a12 * b2;
real d2 = a11 * b2 - b1 * a21;
ct = Point(d1 / d, d2 / d);
r = distance(a, ct);
}
// 面积
real area() const {
return PI * r*r;
}
// 周长
real circumference() const {
return 2 * PI*r;
}
};
// 点与圆的位置关系
// @return -1:点在圆内 0:点在圆上 1:点在圆外
int position(Point p, Circle2D c) {
return dcmp((real)(p - c.ct), c.r);
}
// 直线与圆的位置关系
// @return -1:相交 0:相切 1:相离
int position(Line l, Circle2D c) {
return dcmp(distance(c.ct, l), c.r);
}
// 两圆的位置关系
// @return 0:内含 1:内切 2:相交 3:外切 4:相离
int position(Circle2D a, Circle2D b) {
real d = distance(a.ct, b.ct);
int cmp1 = dcmp(d, a.r + b.r), cmp2 = dcmp(d, abs(a.r - b.r));
if (cmp1 >= 0) return cmp1 + 3;
else return cmp2 + 1;
}
// 圆与直线交点
pair<Point, Point> intersect(Line l, Circle2D c) {
real x0 = c.ct.x, y0 = c.ct.y, r = c.r;
real x1 = l.a.x, y1 = l.a.y;
real x2 = l.b.x, y2 = l.b.y;
real dx = x2 - x1, dy = y2 - y1;
real A = dx * dx + dy * dy;
real B = 2 * dx*(x1 - x0) + 2 * dy*(y1 - y0);
real C = (x1 - x0)*(x1 - x0) + (y1 - y0)*(y1 - y0) - r * r;
real delta = B * B - 4 * A*C;
delta = max((real)0, delta); // 更好地处理相切情况
real t1 = (-B - sqrt(delta)) / 2 / A;
real t2 = (-B + sqrt(delta)) / 2 / A;
return make_pair(Point(x1 + t1 * dx, y1 + t1 * dy), Point(x1 + t2 * dx, y1 + t2 * dy));
}
// @deprecated
// 平面上凸包
// @return 返回上凸包上的点,从左至右,不包含共线点和重合点
vector<Point> hull2DTop(vector<Point> points) {
sort(points.begin(), points.end(), cmpxyz);
vector<Point> rst;
for (const Point& it : points) {
while (rst.size() >= 2) {
if (position(it, Line(rst[rst.size() - 2], rst[rst.size() - 1])) >= 0) {
rst.pop_back();
}
else
break;
}
rst.push_back(it);
}
return rst;
}
// @deprecated
// 平面下凸包
// @return 返回下凸包上的点,从左至右,不包含共线点和重合点
vector<Point> hull2DBottom(vector<Point> points) {
sort(points.begin(), points.end(), cmpxyz);
vector<Point> rst;
for (const Point& it : points) {
while (rst.size() >= 2) {
if (position(it, Line(rst[rst.size() - 2], rst[rst.size() - 1])) <= 0) {
rst.pop_back();
}
else
break;
}
rst.push_back(it);
}
return rst;
}
// 平面凸包
Polygon2D hull2D(vector<Point> points) {
Polygon2D rst(points);
rst.arrange();
return rst;
}

// 平面,由平面上一点和法向量确定
struct Plane {
Point on;
Vector norm;
Plane(Point on = Point(), Vector normal = Vector()) :on(on), norm(normal) {}
// 三点确定平面
Plane(Point a, Point b, Point c) :on(a), norm(cross(b - a, c - a)) {}
};
// 点到平面距离,注意有正负号
real distance(Point p, Plane alpha) {
return (p - alpha.on)*alpha.norm / (real)alpha.norm;
}
// 点与平面位置关系
// @return 0:在平面内 1:在平面正侧 -1:在平面背侧
int position(Point p, Plane alpha) {
return dcmp((p - alpha.on)*alpha.norm);
}
// 直线与平面位置关系
// @return 0:共面 1:平行 2:相交
int position(Line l, Plane alpha) {
if (dcmp(l.direction()*alpha.norm) == 0) {
if (dcmp((l.a - alpha.on)*alpha.norm) == 0)
return 0;
else
return 1;
}
else
return 2;
}
// 点在平面上的投影
Point projection(Point p, Plane alpha) {
return p - zoom(alpha.norm, distance(p, alpha));
}
// 点在直线上的投影
Point projection(Point p, Line l) {
if (colinear(p, l.a, l.b)) return p;
Vector temp = cross(l.direction(), p - l.a);
Point c = l.a + temp;
Plane A(l.a, l.b, c);
return projection(p, A);
}

// 直线左侧表示半平面
typedef Line HalfPlane;

// 半平面交构成凸多边形
Polygon2D intersect(vector<HalfPlane> hps) {
/* 不保证半平面封闭时添加这些
const real INF = 1e50;
hps.emplace_back(Point(-INF, -INF), Point(INF, -INF));
hps.emplace_back(Point(INF, -INF), Point(INF, INF));
hps.emplace_back(Point(INF, INF), Point(-INF, INF));
hps.emplace_back(Point(-INF, INF), Point(-INF, -INF));
*/
// 极角排序
sort(hps.begin(), hps.end(), [](const HalfPlane& a, const HalfPlane& b) {
if (dcmp(a.direction().direction(), b.direction().direction()) != 0)
return a.direction().direction()<b.direction().direction();
else
return position(a.a, b) == 1;
});
// 极角相同保留左侧那一个
auto uend = unique(hps.begin(), hps.end(), [](const HalfPlane& a, const HalfPlane& b) {
return dcmp(a.direction().direction(), b.direction().direction()) == 0;
});
HalfPlane* que = new HalfPlane[uend - hps.begin()];
HalfPlane* frt = que, *bak = que;
for (auto it = hps.begin(); it != uend; ++it) {
while (bak - frt >= 2) {
if (position(intersect(*(bak - 1), *(bak - 2)), *it) == 1)
break;
else
--bak;
}
while (bak - frt >= 2) {
if (position(intersect(*frt, *(frt + 1)), *it) == 1)
break;
else
++frt;
}
*bak++ = *it;
if (bak - frt >= 2 && position(*(bak - 1), *(bak - 2)) != 1) { delete[] que; return Polygon2D(); }
}
// 最后用队首检查队尾
while (bak - frt >= 2) {
if (position(intersect(*(bak - 1), *(bak - 2)), *frt) == 1)
break;
else
--bak;
}
if (bak - frt <= 2) { delete[] que; return Polygon2D(); }
Polygon2D rst;
int n = bak - frt;
rst.vtx.resize(n);
for (int i = 0; i < n; i++) {
rst.vtx[i] = intersect(frt[i], frt[(i + 1) % n]);
}
delete[] que;
return rst;
}
// *凸*多边形交
Polygon2D intersect(Polygon2D a, Polygon2D b) {
int n1 = a.vtx.size(), n2 = b.vtx.size();
vector<HalfPlane> hps(n1 + n2);
for (int i = 0; i < n1; i++) {
hps[i] = HalfPlane(a.vtx[i], a.vtx[(i + 1) % n1]);
}
for (int i = 0; i < n2; i++) {
hps[n1 + i] = HalfPlane(b.vtx[i], b.vtx[(i + 1) % n2]);
}
return intersect(hps);
}
// 圆的切线
pair<Line, Line> tangent(Point p, Circle2D c) {
real sina = c.r / distance(p, c.ct);
sina = min((real)1, sina); // 更好地处理点在圆上的情况
Vector v = c.ct - p;
Vector v1 = rotate(v, asin(sina));
Vector v2 = rotate(v, -asin(sina));
return make_pair(Line(p, p + v1), Line(p, p + v2));
}
// 圆与顶点在圆心的三角形交的面积
real intersect(Point a, Point b, Circle2D c) {
if (position(a, c) <= 0 && position(b, c) <= 0) {
Vector ca = a - c.ct, cb = b - c.ct;
return abs(Cross(ca, cb)) / 2;
}
else if (position(a, c) <= 0 || position(b, c) <= 0) {
if (position(b, c) <= 0) swap(a, b); // a在内部,b在外部
Vector ca = a - c.ct, cb = b - c.ct;
pair<Point, Point> inter = intersect(Line(a, b), c);
Point d = distance(inter.first, b)<distance(inter.second, b) ? inter.first : inter.second; // 取靠近b的交点
return (abs(Cross(ca, d - c.ct)) + angle(cb, d - c.ct)*c.r*c.r) / 2;
}
else {
if (position(Line(a, b), c) < 0) {
pair<Point, Point> inter = intersect(Line(a, b), c);
Point aa = inter.first, bb = inter.second;
if (distance(a, aa)>distance(a, bb)) swap(aa, bb);
Vector ca = a - c.ct, cb = b - c.ct;
Vector caa = aa - c.ct, cbb = bb - c.ct;
return (abs(Cross(caa, cbb)) + (angle(ca, caa) + angle(cb, cbb))*c.r*c.r) / 2;
}
else {
Vector ca = a - c.ct, cb = b - c.ct;
return angle(ca, cb)*c.r*c.r / 2;
}
}
}
}
using namespace Geometry;

读入输出优化

简单版

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
inline char nc() {
static char buf[100000],*p1=buf,*p2=buf;
return p1==p2&&(p2=(p1=buf)+fread(buf,1,sizeof(buf),stdin),p1==p2)?EOF:*p1++;
}
template<typename T> inline void read(T& x) {
char ch;
for(ch=nc();ch<'0'||ch>'9';ch=nc());
for(x=0;ch>='0'&&ch<='9';x=x*10+(ch&0xf),ch=nc());
}
template<typename T> inline bool readEnd(T& x) {
char ch;
for(ch=nc();(ch<'0'||ch>'9')&&ch!=EOF;ch=nc());
if(ch==EOF) return false;
for(x=0;ch>='0'&&ch<='9';x=x*10+(ch&0xf),ch=nc());
return true;
}
template<typename T> inline void readMore(T& x) {
char ch; int sgn=1, k=0;
for(ch=nc();ch<'0'||ch>'9';ch=='-'&&(sgn=-1),ch=nc());
for(x=0;ch>='0'&&ch<='9'||ch=='.';ch=='.'&&(k=1)||(x=x*10+(ch&0xf),k*=10),ch=nc());
x*=sgn; if(k)x/=k;
}

复杂版

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
class IO {
private:
static const int IO_BUFF_SIZE = 1e6;
char buf[IO_BUFF_SIZE], *p, *q;
char outBuf[IO_BUFF_SIZE], *pOut;
inline void flushOut() {
fwrite(outBuf, 1, pOut - outBuf, stdout);
pOut = outBuf;
}
static inline bool blank(char ch) {
return ch == ' ' || ch == '\n' || ch == '\t' || ch == '\r';
}

public:
IO() :p(buf), q(buf), pOut(outBuf) {}
~IO() { flushOut(); }
inline bool next(char& ch) {
if (p == q) {
p = buf;
q = buf + fread(buf, 1, IO_BUFF_SIZE, stdin);
if (q == buf) return false;
}
ch = *p++;
return true;
}
inline bool next(char* str) {
char ch = ' ';
while (next(ch) && blank(ch));
if (blank(ch)) return false;
for (*str++ = ch; next(*str) && !blank(*str); ++str);
*str = '\0';
return true;
}
template<typename T>
inline bool nextNum(T& x) {
char ch = '\0';
int sgn = 1, k=0;
while (next(ch) && !isdigit(ch)) if(ch=='-') sgn = -1;
if (!isdigit(ch)) return false;
for (x = ch & 0xf; next(ch) && (isdigit(ch) || ch == '.'); ) {
if(ch=='.') k=1;
else {x = x * 10 + (ch & 0xf); k*=10;}
}
if(k) x/=k;
x*=sgn;
return true;
}
template<typename T>
inline bool next(T& x) {
char ch = '\0';
while (next(ch) && !isdigit(ch));
if (!isdigit(ch)) return false;
for (x = ch & 0xf; next(ch) && isdigit(ch); x = x * 10 + (ch & 0xf));
return true;
}
inline void print(char ch) {
*pOut++ = ch;
if (outBuf + IO_BUFF_SIZE == pOut) {
flushOut();
}
}
inline void print(const char* str) {
while (*str) print(*str++);
}
inline void print(char* str) {
while (*str) print(*str++);
}
template<typename T>
inline void print(T x) {
if (x < 0) { print('-'); x = -x; }
if (x >= 10) print(x / 10);
print((char)(x % 10 | 0x30));
}
template<typename T>
inline void println(T x) {
print(x);
print('\n');
}
} io;

KMP

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
#include <cstdio>
#include <cstring>

int ns, nss;
char s[10001], ss[1000001];
int Next[10001], fail[10001];

void makeNext() {
int k = 0;
for (int i = 1; i <= ns; i++) {
Next[i] = k;
fail[i] = s[i]==s[k]?fail[k]:k;
while (k && s[i] != s[k])
k = fail[k];
if (s[i] == s[k])
k++;
}
}

// 可重叠查找
int match() {
int rst = 0;
int j = 0;
for (int i = 0; i < nss; i++) {
while (j && ss[i] != s[j])
j = fail[j];
if (ss[i] == s[j]) {
j++;
if (j == ns) {
rst++;
}
}
}
return rst;
}

// 循环节长度
int repetend() {
int len = ns-Next[ns];
return ns%len?ns:len;
}

int main() {
int T;
scanf("%d", &T);
while (T--) {
scanf("%s%s", s, ss);
ns = strlen(s);
nss = strlen(ss);
makeNext();
printf("%d\n", match());
}
return 0;
}

Hungary

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
/// HDU 2444
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std;

const int MAXN = 210;
int n,m;
struct {
int v, next;
} e[MAXN*MAXN];
int head[MAXN], coe;
int type[MAXN];
int x[MAXN], nx, y[MAXN], ny;
int match[MAXN]; // i 匹配 match[i]
bool inchain[MAXN]; // hungary的vis,交替链是否搜索了y中某个点

void addEdge(int u, int v) {
e[coe].v = v;
e[coe].next = head[u];
head[u] = coe++;
}

bool dfs(int u, int t) {
type[u] = t;
if(t==1) x[nx++] = u;
else y[ny++] = u;
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
if(type[v]==type[u])
return false;
if(!type[v]) {
if(!dfs(v, -t))
return false;
}
}
return true;
}

bool divide() {
memset(type, 0, sizeof type);
nx = ny = 0;
for(int i=1; i<=n; i++) {
if(!type[i]) {
if(!dfs(i, 1))
return false;
}
}
return true;
}

bool hgrdfs(int u) {
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
if(!inchain[v]) {
inchain[v] = true;
if(match[v] == -1 || hgrdfs(match[v])) {
// 找到增广路,替换增广路的两种边
match[v] = u;
match[u] = v;
return true;
}
}
}
return false;
}

int hungary() {
int rst = 0;
memset(match, -1, sizeof match);
for(int i=0; i<nx; i++) {
int u = x[i];
if(match[u]==-1) {
memset(inchain, 0, sizeof inchain);
if(hgrdfs(u))
rst++;
}
}
return rst;
}

int main() {
while(scanf("%d%d", &n, &m)==2) {
memset(head, -1, sizeof head);
coe = 0;
for(int i=0; i<m; i++) {
int a,b;
scanf("%d%d", &a, &b);
addEdge(a,b);
addEdge(b,a);
}
if(!divide()) {
puts("No");
}
else {
printf("%d\n", hungary());
}
}
return 0;
}

逆元

扩展欧几里得

要求a与mod互质。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
// ax+by=gcd(a,b)
LL exgcd(LL a, LL b, LL& x, LL& y) {
if(b==0) {
x=1; y=0;
return a;
}
else {
LL g = exgcd(b, a%b, y, x);
y-=a/b*x;
return g;
}
}
LL inverse(LL a) {
LL x,y;
exgcd(a, MOD, x, y);
return (x+MOD)%MOD;
}

费马小定理

要求mod为质数。

1
fastPow(a, MOD-2);

线性递推

要求mod为质数。

1
2
3
4
5
6
7
LL inv[MAXN];
void makeInv() {
inv[1] = 1;
for(int i=2; i<MAXN; i++) {
inv[i] = inv[MOD%i]*(MOD-MOD/i)%MOD;
}
}

线性求阶乘的逆元

1
2
3
4
5
6
7
8
9
10
11
LL fac[MAXN], invfac[MAXN];
void makeFac() {
fac[0] = 1;
for(int i=1; i<MAXN; i++) {
fac[i] = fac[i-1]*i%MOD;
}
invfac[MAXN-1] = fastPow(fac[MAXN-1], MOD-2); // 使用费马小定理要求MOD是质数
for(int i=MAXN-2; i>=0; i--) {
invfac[i] = invfac[i+1]*(i+1)%MOD;
}
}

线性筛同时求质数、欧拉函数、莫比乌斯函数、莫比乌斯函数前缀和

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
int np;
LL prime[MAXN], phi[MAXN], mu[MAXN], smu[MAXN];
bool notp[MAXN];

void linearSieve() {
mu[1] = smu[1] = phi[1] = 1;
for (int i = 2; i < MAXN; i++) {
if (!notp[i]) {
prime[np++] = i;
mu[i] = -1;
phi[i] = i-1;
}
for (int j = 0; j < np; j++) {
if(i*prime[j]>=MAXN) break;
notp[i*prime[j]] = true;
if (i%prime[j] == 0) {
mu[i*prime[j]] = 0;
phi[i*prime[j]] = phi[i]*prime[j];
break;
}
else {
mu[i*prime[j]] = -mu[i];
phi[i*prime[j]] = phi[i]*(prime[j]-1);
}
}
smu[i] = smu[i-1]+mu[i];
}
}

\(1 \leq x \leq n\), \(1 \leq y \leq m\), \(gcd(x,y)=k\)的个数

1
2
3
4
5
6
7
8
9
10
LL nGcdEqK(LL n, LL m, LL k) {
n/=k; m/=k;
LL last, rst = 0;
for (LL i = 1; i <= n && i <= m; i = last + 1) {
last = min(n/(n/i), m/(m/i));
rst += (n/i)*(m/i)%p*(smu[last]-smu[i-1])%p;
rst %= p;
}
return rst;
}

FWT

快速沃尔什变换,计算\(C_i=\sum_{i=j \bigoplus k}{A_j \times B_k}\),其中“\(\bigoplus\)”为按位与、按位或或按位异或,“\(\times\)”就是普通乘法。

注意要把n补到大于最大下标的2的整次幂,因为8|7=15,所以补到8是不够的,要补到16。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
void FWT(LL a[],int n)
{
for(int d=1;d<n;d<<=1)
for(int m=d<<1,i=0;i<n;i+=m)
for(int j=0;j<d;j++)
{
LL x=a[i+j],y=a[i+j+d];
//xor:a[i+j]=x+y,a[i+j+d]=x-y;
//and:a[i+j]=x+y;
//or:a[i+j+d]=x+y;
}
}
void UFWT(LL a[],int n)
{
for(int d=1;d<n;d<<=1)
for(int m=d<<1,i=0;i<n;i+=m)
for(int j=0;j<d;j++)
{
LL x=a[i+j],y=a[i+j+d];
//xor:a[i+j]=(x+y)/2,a[i+j+d]=(x-y)/2;
//and:a[i+j]=x-y;
//or:a[i+j+d]=y-x;
}
}
void solve(LL a[], LL b[], int n)
{
FWT(a,n);
FWT(b,n);
for(int i=0;i<n;i++) a[i]=a[i]*b[i];
UFWT(a,n);
}

pb_ds库

priority_queue

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
#include <cstdio>
#include <ext/pb_ds/priority_queue.hpp>
using namespace std;
/*
__gnu_pbds::priority_queue<element_type, cmp=std::less<>, heap_tag=pairing_heap_tag>
binary_heap_tag 二叉堆
binomial_heap_tag 二项堆
rc_binomial_heap_tag
pairing_heap_tag 配对堆
thin_heap_tag
五种操作:push、pop、modify、erase、join
pairing_heap_tag:push和joinO(1),其余均摊O(logn)
binary_heap_tag:只支持push和pop,均为均摊O(logn)
binomial_heap_tag:push为均摊O(1),其余为O(logn)
rc_binomial_heap_tag:push为O(1),其余为O(logn)
thin_heap_tag:push为O(1),不支持join,其余为O(logn);但是如果只有increase_key,modify均摊O(1)
不支持不是不能用,而是用起来很慢
经过实践检测得到的结论:
Dijkstra算法中应用pairing_heap_tag,速度与手写数据结构相当。
配对堆在绝大多数情况下优于二项堆
只有push,pop和join操作时,二叉堆速度较快
有modify操作时,可以考虑thin_heap_tag或者配对堆,或手写数据结构。
*/
typedef __gnu_pbds::priority_queue<int> Heap;
Heap que;

int main() {
Heap::point_iterator it = que.push(0);
que.push(1); que.push(2);
que.modify(it, 3);
printf("%d\n", que.top());
que.erase(it);
printf("%d\n", que.top());

Heap temp;
temp.push(10);
que.join(temp); // temp合并到que上,然后temp清空,O(logn)
printf("%d\n", que.top());
que.erase(10);
return 0;
}

tree

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
#include <cstdio>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/tree_policy.hpp>

/*
typename Key , typename Mapped ,
typename Cmp_Fn = std :: less <Key >,
typename Tag = rb_tree_tag ,
template <
typename Const_Node_Iterator ,
typename Node_Iterator ,
typename Cmp_Fn_ , typename Allocator_ >
class Node_Update = null_tree_node_update ,
typename Allocator = std :: allocator <char > >
class tree ;

tree的类型,可以是rb_tree_tag,splay_tree_tag,ov_tree_tag
Node_Update:可以为空,也可以用pb_ds自带的tree_order_statistics_node_update,这样这个tree就会获得两个函数find_by_order和order_of_key
iterator find_by_order(size_type order) 找第order+1小的元素的迭代器,如果order太大会返回end()
size_type order_of_key(const_key_reference r_key) :询问这个tree中有多少比r_key小的元素

begin(),end(),size(),empty(),clear(),find(const Key),lower_bound(const Key),upper_bound(const Key),erase(iterator),erase(const Key),insert(const pair<> ),operator[]

如果想改成set,只需要将第二个参数Mapped改为null_type(在4.4.0及以下版本的编译器中应用null_mapped_type)就可以了。此时迭代器指向的类型会从pair变成Key,和set几乎没有区别。
当然还有一些其他用法,如:
void join(tree &other) 把other中所有元素移动到*this上(值域不能相交,否则抛出异常。
void split(const_key_reference r_key, tree &other) 清空other,然后把*this中所有大于r_key的元素移动到other。
*/

// 自定义node update求区间和
template < class Node_CItr , class Node_Itr , class Cmp_Fn , class _Alloc >
struct my_node_update {
virtual Node_CItr node_begin () const = 0;
virtual Node_CItr node_end () const = 0;
typedef int metadata_type ;
inline void operator ()( Node_Itr it , Node_CItr end_it ){
Node_Itr l = it. get_l_child (), r = it. get_r_child ();
int left = 0, right = 0;
if(l != end_it ) left = l. get_metadata ();
if(r != end_it ) right = r. get_metadata ();
const_cast < metadata_type &>( it. get_metadata ())
= left + right + (* it)-> second ;
}
inline int prefix_sum (int x) {
int ans = 0;
Node_CItr it = node_begin ();
while (it != node_end ()) {
Node_CItr l = it. get_l_child (), r = it. get_r_child ();
if( Cmp_Fn ()(x, (* it)-> first )) it = l;
else {
ans += (* it)-> second ;
if(l != node_end ()) ans += l. get_metadata ();
it = r;
}
}
return ans;
}
inline int interval_sum (int l, int r){
return prefix_sum (r) - prefix_sum (l - 1);
}
};
int main() {
__gnu_pbds::tree <int , int , std :: less <int >, __gnu_pbds::rb_tree_tag , my_node_update > T;
T [2] = 100; T [3] = 1000; T [4] = 10000;
printf ("%d\n", T. interval_sum (3, 4));
printf ("%d\n", T. prefix_sum (3));
}

hash_table

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
#include <cstdio>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/hash_policy.hpp>

/*
__gnu_pbds::cc_hash_table<Key, Mapped> 拉链法
__gnu_pbds::gp_hash_table<Key, Mapped> 查探法
*/

__gnu_pbds::gp_hash_table<int, int> mp;

int main() {
mp[1] = 3;
mp[5] = 8;
printf("%d\n", mp[1]);
for(auto it=mp.begin(); it!=mp.end(); ++it) {
printf("%d %d\n", it->first, it->second);
}
return 0;
}

最小费用最大流

最大费用取反。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <queue>
#include <utility>
using namespace std;

const int MAXN = 1e5 + 10, MAXM = 1e6 + 10;
const int INF = 0x3f3f3f3f;

struct Edge {
int v, next, cap, flow, cost;
} e[MAXM];
int head[MAXN], coe;
int pre[MAXN], dis[MAXN];
bool vis[MAXN];
int N;
// 每次初始化,节点编号为0~n-1
void init(int n) {
N = n;
coe = 0;
memset(head, -1, sizeof head);
}
void addEdge(int u, int v, int cap, int cost) {
e[coe].v = v;
e[coe].cap = cap;
e[coe].cost = cost;
e[coe].flow = 0;
e[coe].next = head[u];
head[u] = coe++;
e[coe].v = u;
e[coe].cap = 0;
e[coe].cost = -cost;
e[coe].flow = 0;
e[coe].next = head[v];
head[v] = coe++;
}
bool spfa(int s, int t) {
queue<int> que;
memset(dis, INF, sizeof(int)*N);
memset(vis, 0, sizeof(bool)*N);
memset(pre, -1, sizeof(int)*N);
dis[s] = 0;
vis[s] = true;
que.push(s);
while (!que.empty()) {
int u = que.front();
que.pop();
vis[u] = false;
for (int i = head[u]; i != -1; i = e[i].next) {
int v = e[i].v;
if (e[i].cap > e[i].flow && dis[v] > dis[u] + e[i].cost) {
dis[v] = dis[u] + e[i].cost;
pre[v] = i;
if (!vis[v]) {
vis[v] = true;
que.push(v);
}
}
}
}
if (pre[t] == -1) return false;
else return true;
}
// first: 最大流 second: 最小费用
pair<int, int> minCostMaxflow(int s, int t) {
pair<int, int> ans;
while (spfa(s, t)) {
int Min = INF;
for (int i = pre[t]; i != -1; i = pre[e[i ^ 1].v]) {
if (Min > e[i].cap - e[i].flow)
Min = e[i].cap - e[i].flow;
}
for (int i = pre[t]; i != -1; i = pre[e[i ^ 1].v]) {
e[i].flow += Min;
e[i ^ 1].flow -= Min;
ans.second += e[i].cost * Min;
}
ans.first += Min;
}
return ans;
}

int dep[MAXN];

bool dinic_bfs(int source, int dest) {
memset(dep, -1, sizeof(int)*N);
dep[source] = 0;
queue<int> que;
que.push(source);
while (!que.empty()) {
int i = que.front();
que.pop();
for (int j = head[i]; j != -1; j = e[j].next) {
if (dep[e[j].v] < 0 && e[j].cap > e[j].flow) {
dep[e[j].v] = dep[i] + 1;
que.push(e[j].v);
}
}
}
return dep[dest] > 0;
}

inline int dinic_find(int x, int low, int source, int dest) {
if (low <= 0) return false;
if (x == dest) return low;
int cost = 0;
for (int i = head[x]; i != -1; i = e[i].next) {
if (e[i].cap > e[i].flow && dep[e[i].v] == dep[x] + 1) {
int a = dinic_find(e[i].v, min(low - cost, e[i].cap-e[i].flow), source, dest);
if (a > 0) {
cost += a;
e[i].flow += a;
e[i ^ 1].flow -= a;
if (cost >= low)
break;
}
else {
dep[e[i].v] = -1;
}
}
}
return cost;
}

int dinic(int source, int dest) {
int ans = 0;
while (dinic_bfs(source, dest)) {
int tans;
while (tans = dinic_find(source, INF, source, dest))
ans += tans;
}
return ans;
}

k维最远曼哈顿距离

枚举加减法,代码是两个点集之间的最远距离,一个点集只需要一个Min和Max。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
/// HDU 6435
#include <cstdio>
#include <algorithm>
using namespace std;
typedef long long LL;
const int MAXN = 1e5+10;
const LL INF = 0x3f3f3f3f3f3f3f3f;
int n,m,k;
LL a[MAXN][10], b[MAXN][10];

int main() {
int T;
scanf("%d", &T);
while(T--) {
scanf("%d%d%d", &n, &m, &k);
k++;
for(int i=0; i<n; i++) {
for(int j=0; j<k; j++) {
scanf("%lld", a[i]+j);
}
}
for(int i=0; i<m; i++) {
for(int j=0; j<k; j++) {
scanf("%lld", b[i]+j);
}
b[i][0] = -b[i][0];
}
LL ans = 0;
for(int s=0; s<(1<<k); s++) {
LL Mina=INF, Maxa=-INF, Minb=INF, Maxb=-INF;
for(int i=0; i<n; i++) {
LL t = 0;
for(int j=0; j<k; j++) {
if((1<<j)&s)
t+=a[i][j];
else
t-=a[i][j];
}
Mina = min(Mina, t);
Maxa = max(Maxa, t);
}
for(int i=0; i<m; i++) {
LL t = 0;
for(int j=0; j<k; j++) {
if((1<<j)&s)
t+=b[i][j];
else
t-=b[i][j];
}
Minb = min(Minb, t);
Maxb = max(Maxb, t);
}
#define max3(a,b,c) max(max(a,b),c)
ans = max3(ans, Maxa-Minb, Maxb-Mina);
}
printf("%lld\n", ans);
}
return 0;
}

Java常用库

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
import java.io.*;
import java.math.*;
import java.util.*;

public class Main {
static Scanner in = new Scanner(System.in);
public static void main(String[] args) {
int a = in.nextInt();
BigInteger b = in.nextBigInteger();
BigDecimal c = in.nextBigDecimal();
/*
BigDecimal:
构造方法:
BigDecimal(BigInteger val)
BigDecimal(BigInteger unscaledVal, int scale)
BigDecimal(BigInteger unscaledVal, int scale, MathContext mc)
BigDecimal(BigInteger val, MathContext mc)
BigDecimal(char[] in)
BigDecimal(char[] in, int offset, int len)
BigDecimal(char[] in, int offset, int len, MathContext mc)
BigDecimal(char[] in, MathContext mc)
BigDecimal(double val)
BigDecimal(double val, MathContext mc)
BigDecimal(int val)
BigDecimal(int val, MathContext mc)
BigDecimal(long val)
BigDecimal(long val, MathContext mc)
BigDecimal(String val)
BigDecimal(String val, MathContext mc)
成员方法:
BigDecimal abs()
BigDecimal abs(MathContext mc)
BigDecimal add(BigDecimal augend)
BigDecimal add(BigDecimal augend, MathContext mc)
byte byteValueExact()
int compareTo(BigDecimal val)
BigDecimal divide(BigDecimal divisor)
BigDecimal divide(BigDecimal divisor, int roundingMode)
BigDecimal divide(BigDecimal divisor, int scale, int roundingMode)
BigDecimal divide(BigDecimal divisor, int scale, RoundingMode roundingMode)
BigDecimal divide(BigDecimal divisor, MathContext mc)
BigDecimal divide(BigDecimal divisor, RoundingMode roundingMode)
BigDecimal[] divideAndRemainder(BigDecimal divisor)
BigDecimal[] divideAndRemainder(BigDecimal divisor, MathContext mc)
BigDecimal divideToIntegralValue(BigDecimal divisor)
BigDecimal divideToIntegralValue(BigDecimal divisor, MathContext mc)
double doubleValue()
boolean equals(Object x)
float floatValue()
int hashCode()
int intValue()
int intValueExact()
long longValue()
long longValueExact()
BigDecimal max(BigDecimal val)
BigDecimal min(BigDecimal val)
BigDecimal movePointLeft(int n)
BigDecimal movePointRight(int n)
BigDecimal multiply(BigDecimal multiplicand)
BigDecimal multiply(BigDecimal multiplicand, MathContext mc)
BigDecimal negate()
BigDecimal negate(MathContext mc)
BigDecimal plus()
BigDecimal plus(MathContext mc)
BigDecimal pow(int n)
BigDecimal pow(int n, MathContext mc)
int precision()
BigDecimal remainder(BigDecimal divisor)
BigDecimal remainder(BigDecimal divisor, MathContext mc)
BigDecimal round(MathContext mc)
int scale()
BigDecimal scaleByPowerOfTen(int n)
BigDecimal setScale(int newScale)
Returns a BigDecimal whose scale is the specified value, and whose value is numerically equal to this BigDecimal's.
BigDecimal setScale(int newScale, int roundingMode)
BigDecimal setScale(int newScale, RoundingMode roundingMode)
Returns a BigDecimal whose scale is the specified value, and whose unscaled value is determined by multiplying or dividing this BigDecimal's unscaled value by the appropriate power of ten to maintain its overall value.
short shortValueExact()
int signum()
Returns the signum function of this BigDecimal. (1,0,-1)
BigDecimal stripTrailingZeros()
Returns a BigDecimal which is numerically equal to this one but with any trailing zeros removed from the representation.
BigDecimal subtract(BigDecimal subtrahend)
BigDecimal subtract(BigDecimal subtrahend, MathContext mc)
BigInteger toBigInteger()
BigInteger toBigIntegerExact()
String toEngineeringString()
String toPlainString()
String toString()
BigDecimal ulp()
Returns the size of an ulp, a unit in the last place, of this BigDecimal.
BigInteger unscaledValue()
static BigDecimal valueOf(double val)
static BigDecimal valueOf(long val)
static BigDecimal valueOf(long unscaledVal, int scale)
*/
BigDecimal test = new BigDecimal("1.234567");
test = test.setScale(3, RoundingMode.HALF_UP);
System.out.println(test);
test = test.setScale(7, RoundingMode.HALF_EVEN);
System.out.println(test);
test = test.divide(new BigDecimal("3"), MathContext.UNLIMITED); // 默认也是UNLIMITED精度,无限小数会报错
System.out.println(test);

/*
BigInteger:
构造方法:
BigInteger(byte[] val)
BigInteger(int signum, byte[] magnitude)
BigInteger(int bitLength, int certainty, Random rnd)
Constructs a randomly generated positive BigInteger that is probably prime, with the specified bitLength.
BigInteger(int numBits, Random rnd)
Constructs a randomly generated BigInteger, uniformly distributed over the range 0 to (2numBits - 1), inclusive.
BigInteger(String val)
BigInteger(String val, int radix)
成员方法:
BigInteger abs()
BigInteger add(BigInteger val)
BigInteger and(BigInteger val)
BigInteger andNot(BigInteger val)
Returns a BigInteger whose value is (this & ~val).
int bitCount()
Returns the number of bits in the two's complement representation of this BigInteger that differ from its sign bit.
int bitLength()
Returns the number of bits in the minimal two's-complement representation of this BigInteger, excluding a sign bit.
BigInteger clearBit(int n)
Returns a BigInteger whose value is equivalent to this BigInteger with the designated bit cleared.
int compareTo(BigInteger val)
BigInteger divide(BigInteger val)
BigInteger[] divideAndRemainder(BigInteger val)
double doubleValue()
boolean equals(Object x)
BigInteger flipBit(int n)
Returns a BigInteger whose value is equivalent to this BigInteger with the designated bit flipped.
float floatValue()
BigInteger gcd(BigInteger val)
Returns a BigInteger whose value is the greatest common divisor of abs(this) and abs(val).
int getLowestSetBit()
Returns the index of the rightmost (lowest-order) one bit in this BigInteger (the number of zero bits to the right of the rightmost one bit).
int hashCode()
int intValue()
boolean isProbablePrime(int certainty)
Returns true if this BigInteger is probably prime, false if it's definitely composite.
long longValue()
BigInteger max(BigInteger val)
BigInteger min(BigInteger val)
BigInteger mod(BigInteger m)
BigInteger modInverse(BigInteger m)
Returns a BigInteger whose value is (this^-1 mod m).
BigInteger modPow(BigInteger exponent, BigInteger m)
BigInteger multiply(BigInteger val)
BigInteger negate()
BigInteger nextProbablePrime()
Returns the first integer greater than this BigInteger that is probably prime.
BigInteger not()
BigInteger or(BigInteger val)
BigInteger pow(int exponent)
static BigInteger probablePrime(int bitLength, Random rnd)
Returns a positive BigInteger that is probably prime, with the specified bitLength.
BigInteger remainder(BigInteger val)
Returns a BigInteger whose value is (this % val).
BigInteger setBit(int n)
Returns a BigInteger whose value is equivalent to this BigInteger with the designated bit set.
BigInteger shiftLeft(int n)
Returns a BigInteger whose value is (this << n).
BigInteger shiftRight(int n)
Returns a BigInteger whose value is (this >> n).
int signum()
BigInteger subtract(BigInteger val)
boolean testBit(int n)
Returns true if and only if the designated bit is set.
byte[] toByteArray()
String toString()
String toString(int radix)
static BigInteger valueOf(long val)
BigInteger xor(BigInteger val)
*/

MyPair[] pairs = new MyPair[1000];
Arrays.sort(pairs);
Arrays.binarySearch(pairs, 1, 4, new MyPair());
/*
二分查找
如果元素在数组中,则值为0~n-1,否则值为-1~-(n+1),表示第一个比它大的值的位置,下标从1开始
*/
List<MyPair> pairList = new ArrayList<>();
pairList.add(new MyPair());
//pairList.sort();
pairList.sort(new Cmp());
Collections.shuffle(pairList);
Collections.swap(pairList, 1, 3);
Collections.sort(pairList);
}
}

class MyPair implements Comparable {
int x, y;

@Override
public int compareTo(Object o) {
MyPair b = (MyPair)o;
if(x!=b.x) return x<b.x?-1:1;
else if(y!=b.y) return y<b.y?-1:1;
return 0;
}
}

class Cmp implements Comparator<MyPair> {
@Override
public int compare(MyPair o1, MyPair o2) {
if(o1.x!=o2.x) return o1.x<o2.x?-1:1;
else if(o1.y!=o2.y) return o1.y<o2.y?-1:1;
return 0;
}
}

对拍

check.bat (Windows BAT)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
@echo off
echo building rand...
g++ rand.cpp -O2 -std=c++14 -Wall -o rand.exe
echo building std...
g++ std.cpp -O2 -std=c++14 -Wall -o std.exe
echo building mine...
g++ mine.cpp -O2 -std=c++14 -Wall -o mine.exe
echo start checking...
:loop
rand.exe %random% > in.txt
std.exe < in.txt > stdout.txt
mine.exe < in.txt > mineout.txt
fc stdout.txt mineout.txt
if not errorlevel 1 goto loop
pause
goto loop

check.sh (Linux shell)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#!/bin/bash
set -e
echo building rand...
g++ rand.cpp -std=c++14 -O2 -Wall -o rand.o
echo building std...
g++ std.cpp -std=c++14 -O2 -Wall -o std.o
echo building mine...
g++ mine.cpp -std=c++14 -O2 -Wall -o mine.o
echo start checking...
while true; do
./rand.o $RANDOM > in.txt
./mine.o < in.txt > mineout.txt
./std.o < in.txt > stdout.txt
diff stdout.txt mineout.txt
echo OK
done

rand.cpp

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cstdlib>
#include <random>
using namespace std;

default_random_engine e;

inline int randInt(int min, int max) {
return uniform_int_distribution<int>(min, max)(e);
}

int main(int argc, char const *argv[]) {
e.seed(strtoul(argv[1], nullptr, 0));

return 0;
}
/*
#include <cstdio>
#include <random>
#include <functional>
#include <chrono>
#include <cstdlib>
using namespace std;

long long seed1 = chrono::nanoseconds( chrono::system_clock::now().time_since_epoch() ).count();
default_random_engine e(seed1);
uniform_int_distribution<int> dis10(1, 10);

inline int randInt(int min, int max) {
unsigned temp = (rand() << 15 | rand()) << 2 | (rand() & 3);
return temp%(unsigned)(max-min+1)+min;
}

int main(int argc, char const *argv[]) {
srand(strtol(argv[1], nullptr, 0));

for(int i=0; i<10; i++) {
printf("%d\n", randInt(-1e9, 1e9));
}
return 0;
}
*/

FFT

快速傅里叶变换,计算\(C_i=\sum_{i=j + k}{A_j \times B_k}\),即普通的多项式相乘。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <complex>
using namespace std;

typedef complex<double> Complex;
const double PI = acos(-1);
const int MAXN = (1<<20)+10;
int N; // 大于等于结果项数(n1+n2-1)的2的整次幂,init会赋值
Complex omega[MAXN], omegaInv[MAXN];
int n1, n2;
double a[MAXN], b[MAXN], c[MAXN];

void init(int n1, int n2) {
for(N=1; N<(n1+n2-1); N<<=1);
for(int i=0; i<N; i++) {
omega[i] = Complex(cos(2*PI*i/N), sin(2*PI*i/N));
omegaInv[i] = conj(omega[i]);
}
}

void dft(Complex a[], Complex ome[] = omega) {
for(int i=0,k=log2(N); i<N; i++) {
int t = 0;
for(int j=0; j<k; j++) {
if(i&(1<<j))
t|=1<<(k-j-1);
}
if(i<t) swap(a[i], a[t]);
}
for(int l=2; l<=N; l<<=1) {
int m = l>>1;
for(Complex* p=a; p<a+N; p+=l) {
for(int i=0; i<m; i++) {
Complex t = ome[N/l*i]*p[m+i];
p[m+i]=p[i]-t;
p[i]+=t;
}
}
}
}

void idft(Complex a[]) {
dft(a, omegaInv);
for(int i=0; i<N; i++)
a[i]/=N;
}

// 多项式相乘,f(x)=a[i]*x^i, g(x)=b[i]*x^i, i=0..n-1,n1, n2是项数
void multiply(const double a[], int n1, const double b[], int n2, double rst[]) {
static Complex ia[MAXN], ib[MAXN];
init(n1, n2);
copy(a, a+n1, ia);
fill(ia+n1, ia+N, Complex());
copy(b, b+n2, ib);
fill(ib+n2, ib+N, Complex());
dft(ia); dft(ib);
for(int i=0; i<N; i++)
ia[i]*=ib[i];
idft(ia);
for(int i=0; i<n1+n2-1; i++)
rst[i] = ia[i].real();
}

int main() {
while(scanf("%d%d", &n1, &n2)==2) {
for(int i=0; i<n1; i++) {
scanf("%lf", a+i);
}
for(int i=0; i<n2; i++) {
scanf("%lf", b+i);
}
multiply(a, n1, b, n2, c);
for(int i=0; i<n1+n2-1; i++) {
printf("%f%c", c[i], " \n"[i==n1+n2-2]);
}
}
return 0;
}

数值积分

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
const double EPS = 1e-10;

double F(double x) {
return x*x*x+2*x*x+0.4*x+3.332;
}

// 阶数为4或5用科特斯公式
double cotes(double a, double b) {
double d=(a+b)/2, c=(a+d)/2, e=(d+b)/2;
return (7*F(a)+32*F(c)+12*F(d)+32*F(e)+7*F(b))*(b-a)/90;
}
// 阶数小于等于3用辛普森公式
double simpson(double a, double b) {
double c = (a+b)/2;
return (F(a)+4*F(c)+F(b))*(b-a)/6;
}

double asr(double a, double b, double eps, double A) {
double c = (a+b)/2;
double L = simpson(a, c), R = simpson(c, b);
//double L = cotes(a, c), R = cotes(c, b);
if(abs(L+R-A) <= 15*eps) return L+R+(L+R-A)/15;
return asr(a, c, eps/2, L) + asr(c, b, eps/2, R);
}
// 调用这个
double asr(double a, double b, double eps){
return asr(a, b, eps, simpson(a,b));
//return asr(a, b, eps, cotes(a,b));
}

Dijkstra+pb_ds

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
/// 洛谷P4779
#include <cstdio>
#include <cstring>
#include <ext/pb_ds/priority_queue.hpp>
using namespace std;
typedef long long LL;
const int MAXN = 1e5+10, MAXM = 1e5+10;
const LL INF = 0x3f3f3f3f3f3f3f3f;
int n, m;
struct {
int v,next;
LL w;
} e[MAXM*2];
int head[MAXN], coe;
LL mincost[MAXN]; // dijkstra result

struct Dis {
int u;
LL d;
bool operator< (const Dis& b) const {
return b.d<d;
}
};
typedef __gnu_pbds::priority_queue<Dis> Heap;

inline void addEdge(int u, int v, LL w) {
e[coe] = {v,head[u],w};
head[u] = coe++;
}

void dijkstra(int source, int dest = -1) {
static Heap::point_iterator hit[MAXN];
//static bool vis[MAXN];
static bool hed[MAXN];

Heap que;
memset(mincost, 0x3f, sizeof mincost);
//memset(vis, 0, sizeof vis);
memset(hed, 0, sizeof hed);
mincost[source] = 0;
hit[source] = que.push({source, 0});
hed[source] = true;

while(!que.empty()) {
int u = que.top().u;
LL d = que.top().d;
if(u==dest) return;
que.pop();
hed[u] = false;
//vis[u] = true;
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
LL w = e[i].w;
//if(vis[v]) continue;
if(d+w<mincost[v]) {
mincost[v] = d+w;
if(hed[v])
que.modify(hit[v], {v, d+w});
else {
hit[v] = que.push({v, d+w});
hed[v] = true;
}
}
}
}
}

int main() {
int s;
while(scanf("%d%d%d", &n, &m, &s)==3) {
memset(head, -1, sizeof head);
coe = 0;
for(int i=0; i<m; i++) {
int a,b,c;
scanf("%d%d%d", &a, &b, &c);
addEdge(a,b,c);
}
dijkstra(s);
for(int i=1; i<=n; i++) {
printf("%lld%c", mincost[i], " \n"[i==n]);
}
}
return 0;
}

AC自动机

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
/// HDU 3065
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <queue>
using namespace std;

const int MAXN = 50 * 1000 + 10;
int n;
char mode[1010][60];
char a[2000000 + 10];
int ans[1010];

struct Node {
int val;
int fail;
int next[26];
} node[MAXN];
int con;

void insertTrie(char str[], int id) {
int p = 0;
for (int i = 0; str[i]; i++) {
int ha = str[i] - 'A';
if (!node[p].next[ha]) {
node[p].next[ha] = con++;
}
p = node[p].next[ha];
}
node[p].val = id;
}

void makeFail() {
queue<int> que;
que.push(0);
while (!que.empty()) {
int p = que.front();
que.pop();
for (int i = 0; i<26; i++) {
if (!node[p].next[i]) continue;
if (p == 0) node[node[p].next[i]].fail = 0;
else {
int temp = node[p].fail;
while (temp) {
if (node[temp].next[i]) {
node[node[p].next[i]].fail = node[temp].next[i];
break;
}
temp = node[temp].fail;
}
if (!temp)
node[node[p].next[i]].fail = node[0].next[i];
}
que.push(node[p].next[i]);
}
}
}

void match(char str[]) {
int p = 0;
for (int i = 0; str[i]; i++) {
int ha = str[i] - 'A';
while (p && (!(ha >= 0 && ha<26) || !node[p].next[ha]))
p = node[p].fail;
if (ha >= 0 && ha<26 && node[p].next[ha])
p = node[p].next[ha];
int temp = p;
while (temp) {
if (node[temp].val)
ans[node[temp].val]++;
temp = node[temp].fail;
}
}
}

void print(int p) {
printf("%d:\n", p);
printf("val: %d\n", node[p].val);
printf("fail: %d\n", node[p].fail);
printf("edge:\n");
for (int i = 0; i < 26; i++) {
if (node[p].next[i]) {
printf("%c: %d\n", 'A'+i, node[p].next[i]);
}
}
printf("\n");
for (int i = 0; i < 26; i++) {
if (node[p].next[i]) {
print(node[p].next[i]);
}
}
}

int main() {
while (scanf("%d", &n) == 1) {
memset(ans, 0, sizeof ans);
memset(node, 0, sizeof node);
con = 1;
for (int i = 1; i <= n; i++) {
scanf("%s", mode[i]);
insertTrie(mode[i], i);
}
makeFail();
//print(0);
scanf("%s", a);
match(a);
for (int i = 1; i <= n; i++) {
if (ans[i]) {
printf("%s: %d\n", mode[i], ans[i]);
}
}
}
return 0;
}

球盒问题

n个球 m个盒 空盒 情况数(dp、S的含义见下方)
相同 相同 允许 dp[n][m]
相同 相同 不允许 n>=m时dp[n-m][m],否则为0
相同 不同 允许 \(\binom{n+m-1}{m-1}\)
相同 不同 不允许 \(\binom{n-1}{m-1}\)
不同 相同 允许 \(\sum_{i=0}^{m}{S(n,i)}\)
不同 相同 不允许 \(S(n,m)\)
不同 不同 允许 \(m^n\)
不同 不同 不允许 \(m! \times S(n,m)\)

第二类斯特林数

OEIS A008277

S[i][j]的含义是,把大小为i的集合划分为j个非空子集合的方案数。

1
2
3
S[i][i]=1, i>=0 // i个球放入i个盒子
S[i][j]=S[i-1][j]*j+S[i-1][j-1], i>=2, 1<=j<=i-1 // 前i-1个球在j个盒子中,第i个球随意放;或前i-1个球在j-1个盒子中,第i个球只能放在第j个盒子中
其余情况为0
1
2
3
4
5
6
7
8
9
10
11
LL S[MAXN][MAXN];
void init() {
for(int i=0; i<MAXN; i++) {
S[i][i] = 1;
}
for(int i=2; i<MAXN; i++) {
for(int j=1; j<i; j++) {
S[i][j]=S[i-1][j]*j+S[i-1][j-1];
}
}
}
1
2
3
4
5
6
7
8
9
10
11
1
1 1
1 3 1
1 7 6 1
1 15 25 10 1
1 31 90 65 15 1
1 63 301 350 140 21 1
1 127 966 1701 1050 266 28 1
1 255 3025 7770 6951 2646 462 36 1
1 511 9330 34105 42525 22827 5880 750 45 1
1 1023 28501 145750 246730 179487 63987 11880 1155 55 1

球、盒均相同的DP方法

OEIS A026820

dp[i][j]的含义是,把非负整数i划分为j个非负整数之和的方案数。

1
2
3
4
5
6
dp[0][i]=1, i>=0 // 0个球放在i个盒子中
dp[1][i]=1, i>=1 // 1个球放在i个盒子中
dp[i][1]=1, i>=0 // i个球放在1个盒子中
dp[i][j]=dp[i][j-1]+dp[i-j][j], i>=2, j<=i // i个球放在j-1个盒子中;或i-j个球放在j个盒子中,再给每个盒子里放一个球
dp[i][j]=dp[i][j-1], i>=2, j>i // i个球放在j-1个盒子中
其余情况为0
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
LL dp[MAXN][MAXN];
void init() {
for(int i=0; i<MAXN; i++)
dp[0][i] = 1;
for(int i=1; i<MAXN; i++)
dp[1][i] = 1;
for(int i=0; i<MAXN; i++)
dp[i][1] = 1;
for(int i=2; i<MAXN; i++) {
for(int j=2; j<MAXN; j++) {
dp[i][j] = dp[i][j-1];
if(i>j) dp[i][j] += dp[i-j][j];
}
}
}
1
2
3
4
5
6
7
8
9
10
11
12
1
1 2
1 2 3
1 3 4 5
1 3 5 6 7
1 4 7 9 10 11
1 4 8 11 13 14 15
1 5 10 15 18 20 21 22
1 5 12 18 23 26 28 29 30
1 6 14 23 30 35 38 40 41 42
1 6 16 27 37 44 49 52 54 55 56
1 7 19 34 47 58 65 70 73 75 76 77

第一类斯特林数

S[i][j]的含义是,把i个物品排成j个非空循环排列的方案数。

1
2
3
S[0][0]=1
S[i][j]=(i-1)*S[i-1][j]+S[i-1][j-1], i>=1, j>=1
其余为0
1
2
3
4
5
6
7
8
9
LL S[MAXN][MAXN];
void init() {
S[0][0] = 1;
for(int i=1; i<MAXN; i++) {
for(int j=1; j<=i; j++) {
S[i][j] = (i-1)*S[i-1][j]+S[i-1][j-1];
}
}
}
1
2
3
4
5
6
7
8
9
10
1
1 1
2 3 1
6 11 6 1
24 50 35 10 1
120 274 225 85 15 1
720 1764 1624 735 175 21 1
5040 13068 13132 6769 1960 322 28 1
40320 109584 118124 67284 22449 4536 546 36 1
362880 1026576 1172700 723680 269325 63273 9450 870 45 1

Catalan数

OEIS A000108 \[ C_n=\frac{1}{n+1}\binom{2n}{n}=\frac{(2n)!}{(n+1)!n!}=\binom{2n}{n}-\binom{2n}{n+1} \]

1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862, 16796, 58786, 208012, 742900, 2674440, 9694845, 35357670, 129644790, 477638700, 1767263190

组合数学中有非常多的组合结构可以用卡塔兰数来计数。在Richard P. Stanley的Enumerative Combinatorics: Volume 2一书的习题中包括了66个相异的可由卡塔兰数表达的组合结构。

  • \(C_n\)表示长度2n的dyck word的个数。Dyck word是一个有n个X和n个Y组成的字串,且所有的前缀字串皆满足X的个数大于等于Y的个数。以下为长度为6的dyck words:
XXXYYY XYXXYY XYXYXY XXYYXY XXYXYY
  • 将上例的X换成左括号,Y换成右括号,\(C_n\)表示所有包含n组括号的合法运算式的个数:
((())) ()(()) ()()() (())() (()())
  • \(C_n\)表示有n个节点组成不同构二叉树的方案数。下图中,n等于3,圆形表示节点,月牙形表示什么都没有。

  • \(C_n\)表示有2n+1个节点组成不同构满二叉树的方案数。下图中,n等于3,圆形表示内部节点,月牙形表示外部节点。本质同上。

Catalan number binary tree example.png

  • \(C_n\)表示所有在n×n格点中不越过对角线的单调路径的个数。一个单调路径从格点左下角出发,在格点右上角结束,每一步均为向上或向右。计算这种路径的个数等价于计算Dyck word的个数:X代表“向右”,Y代表“向上”。下图为n=4的情况:

Catalan number 4x4 grid example.svg

  • \(C_n\)表示通过连结顶点而将n+2边的凸多边形的方法个数。下图中为n=4的情况:

Catalan-Hexagons-example.svg

  • \(C_n\)表示对{1, ..., n}依序进出栈的置换个数。一个置换w是依序进出栈的当S(w)=(1, ..., n),其中S(w)递归定义如下:令w=unv,其中n为w的最大元素,u和v为更短的数列;再令S(w)=S(u)S(v)n,其中S为所有含一个元素的数列的单位元。

  • \(C_n\)表示集合{1, ..., n}的不交叉划分的个数.那么, \(C_n\)永远不大于第n项贝尔数. \(C_n\)也表示集合{1, ..., 2n}的不交叉划分的个数,其中每个段落的长度为2。综合这两个结论,可以用数学归纳法证明:在魏格纳半圆分布定律中度数大于2的情形下,所有自由的累积量s为零。 该定律在自由概率论和随机矩阵理论中非常重要。

  • \(C_n\)表示用n个长方形填充一个高度为n的阶梯状图形的方法个数。下图为n=4的情况:

Catalan stairsteps 4.svg

  • \(C_n\)表示表为2×n的矩阵的标准杨氏矩阵的数量。也就是说,它是数字 1, 2, ..., 2n被放置在一个2×n的矩形中并保证每行每列的数字升序排列的方案数。同样的,该式可由勾长公式的一个特殊情形推导得出。

  • \(C_n\)表示n个无标号物品的半序的个数。

错排

考虑一个有n个元素的排列,若一个排列中所有的元素都不在自己原来的位置上,那么这样的排列就称为原排列的一个错排。 n个元素的错排数记为D(n)。 \[ \begin{aligned} D(n)&=\left\{ \begin{aligned} & 0 & n=1 \\ & 1 & n=2 \\ & (n-1)[D(n-1)+D(n-2)] & n \geq 3 \end{aligned} \right. \\ &=n![\frac1{2!}-\frac1{3!}+\dots+(-1)^n\frac1{n!}] \\ &=\lfloor\frac{n!}e+0.5\rfloor \end{aligned} \]

康拓展开

X表示一个排列在所有的全排列中排第几个(从0开始)。 \[ X=A_n(n−1)!+A_{n−1}(n−2)!+\dots+A_1\times0! \] 其中\(A_i\)表示这个排列里从左到右第i个数字之后有多少比这个数字小的数字。

逆运算

假设求4位数中第18个位置(第0个是1 2 3 4)的数字。

18对3!作除法→得3余0

0对2!作除法→得0余0

0对1!作除法→得0余0

据上面的可知:

我们第一位数(最左面的数),比第一位数小的数有3个,显然第一位数为4。

比第二位数小的数字有0个,所以第二位数为1。

比第三位数小的数字有0个,因为1已经用过,所以第三位数为2。

第四位数剩下3。

该数字为 4123。

rope库

1
2
#include <ext/rope>
typedef __gnu_cxx::rope<int> List;
函数 功能
push_back(x) 在末尾添加x,也支持类似的其他deque操作
replace(pos, len, begin, end) 将pos~pos+len-1的值替换为begin~end-1
replace(pos, len, begin, count) 将pos~pos+len-1的值替换为begin~begin+count-1
insert(pos, begin, end) 在pos之前插入,也有其他类似replace的用法
copy(pos, len, begin) replace的等长版
substr(pos, len) 截取
拷贝构造函数 可持久化的方法,O(1)

莫比乌斯反演

F(n)和f(n)是定义在非负整数集合上的两个函数。

定理1

如果: \[ F(n)=\sum_{d|n(d是n的因数)}{f(d)} \] 那么: \[ f(n)=\sum_{d|n}{\mu(d)F(\frac nd)} \]

定理2

如果: \[ F(n)=\sum_{n|d(d是n的倍数)}{f(d)} \] 那么: \[ f(n)=\sum_{n|d}{\mu(\frac dn)F(d)} \]

数位DP

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
typedef long long LL;
int a[20];
LL dp[20][state];//不同题目状态不同
LL dfs(int pos,/*state变量*/,bool lead/*前导零*/,bool limit/*数位上界变量*/)//不是每个题都要判断前导零
{
//递归边界,既然是按位枚举,最低位是0,那么pos==-1说明这个数我枚举完了
if(pos==-1) return 1;/*这里一般返回1,表示你枚举的这个数是合法的,那么这里就需要你在枚举时必须每一位都要满足题目条件,也就是说当前枚举到pos位,一定要保证前面已经枚举的数位是合法的。不过具体题目不同或者写法不同的话不一定要返回1 */
//第二个就是记忆化(在此前可能不同题目还能有一些剪枝)
if(!limit && !lead && dp[pos][state]!=-1) return dp[pos][state];
/*常规写法都是在没有限制的条件记忆化,这里与下面记录状态是对应,具体为什么是有条件的记忆化后面会讲*/
int up=limit?a[pos]:9;//根据limit判断枚举的上界up;这个的例子前面用213讲过了
LL ans=0;
//开始计数
for(int i=0;i<=up;i++)//枚举,然后把不同情况的个数加到ans就可以了
{
if() ...
else if()...
ans+=dfs(pos-1,/*状态转移*/,lead && i==0,limit && i==a[pos]) //最后两个变量传参都是这样写的
/*这里还算比较灵活,不过做几个题就觉得这里也是套路了
大概就是说,我当前数位枚举的数是i,然后根据题目的约束条件分类讨论
去计算不同情况下的个数,还有要根据state变量来保证i的合法性,比如题目
要求数位上不能有62连续出现,那么就是state就是要保存前一位pre,然后分类,
前一位如果是6那么这意味就不能是2,这里一定要保存枚举的这个数是合法*/
}
//计算完,记录状态
if(!limit && !lead) dp[pos][state]=ans;
/*这里对应上面的记忆化,在一定条件下时记录,保证一致性,当然如果约束条件不需要考虑lead,这里就是lead就完全不用考虑了*/
return ans;
}
LL solve(LL x)
{
int pos=0;
while(x)//把数位都分解出来
{
a[pos++]=x%10;//个人老是喜欢编号为[0,pos),看不惯的就按自己习惯来,反正注意数位边界就行
x/=10;
}
return dfs(pos-1/*从最高位开始枚举*/,/*一系列状态 */,true,true);//刚开始最高位都是有限制并且有前导零的,显然比最高位还要高的一位视为0嘛
}
int main()
{
LL le,ri;
while(~scanf("%LLd%LLd",&le,&ri))
{
//初始化dp数组为-1,这里还有更加优美的优化,后面讲
printf("%LLd\n",solve(ri)-solve(le-1));
}
}

大整数类

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
#include <cstdio>
#include <algorithm>
#include <iostream>
#include <string>
#include <cstring>
#include <vector>
#include <cctype>
using namespace std;
typedef long long LL;

struct BigInteger {
static const int D = 10000;
typedef vector<int> DigSet;

enum Sign {positive, negative} sign;
DigSet dig;
BigInteger(const DigSet& dig = DigSet(), Sign sign = positive):sign(sign),dig(dig){}
BigInteger(LL val):sign(val>=0?positive:negative) {
val = abs(val);
while (val) {
dig.push_back(val%D);
val /= D;
}
}
BigInteger(const string& val) {
for (int i = val.size() - 1; i >= 0; i-=4) {
int d = 0;
for (int j = 0, p=1; j < 4 && i - j >= 0 && isdigit(val[i - j]); j++, p*=10) {
d += (val[i-j] & 0xf)*p;
}
dig.push_back(d);
}
sign = Sign(val[0] == '-');
purify();
}
string toString() const {
if (zero()) return "0";
string rst;
if (sign==negative) rst = "-";
rst += to_string(highest());
for (int i = 2; i <= dig.size(); i++) {
static char buf[10];
sprintf(buf, "%04d", highest(i));
rst += buf;
}
return rst;
}
friend istream& operator>> (istream& is, BigInteger& x) {
string s;
is >> s;
x = s;
return is;
}
friend ostream& operator<< (ostream& os, const BigInteger& x) {
os << x.toString();
return os;
}

int highest(int i=1) const {
return dig[dig.size() - i];
}
bool zero() const {
return dig.empty();
}
BigInteger& purify() {
if (zero()) {
sign = positive;
return *this;
}
for (int i = 0; i < dig.size()-1; i++) {
dig[i + 1] += dig[i] / D;
dig[i] %= D;
if (dig[i] < 0) {
dig[i + 1]--;
dig[i] += D;
}
}
while (highest() >= D) {
dig.push_back(highest() / D);
dig[dig.size() - 2] %= D;
}
while (!zero() && !highest()) {
dig.pop_back();
}
if (zero()) {
sign = positive;
}
return *this;
}


static DigSet add(const DigSet& x, const DigSet& y) {
DigSet rst = x.size() > y.size() ? x : y;
const DigSet& less = x.size() > y.size() ? y : x;
for (int i = 0; i < less.size(); i++) {
rst[i] += less[i];
}
return rst;
}
static DigSet sub(DigSet x, const DigSet& y) {
for (int i = 0; i < y.size(); i++) {
x[i] -= y[i];
}
return x;
}
static int cmp(const DigSet& x, const DigSet& y) {
if (x.size() < y.size()) return -1;
else if (x.size() > y.size()) return 1;
for (int i = x.size()-1; i >= 0; i--) {
if (x[i] < y[i]) return -1;
else if (x[i] > y[i]) return 1;
}
return 0;
}
bool operator< (const BigInteger& b) const {
if (sign != b.sign) return sign==negative;
else return cmp(dig, b.dig) == (int)sign*2-1;
}
bool operator== (const BigInteger& b) const {
return sign == b.sign && cmp(dig, b.dig) == 0;
}
bool operator<= (const BigInteger& b) const {
return *this == b || *this < b;
}
BigInteger operator+ (const BigInteger& b) const {
if (sign == b.sign) {
return BigInteger(add(dig, b.dig), sign).purify();
}
else {
return *this- -b;
}
}
BigInteger operator-() const {
return BigInteger(dig, (Sign)!sign);
}
BigInteger operator-(const BigInteger& b) const {
if (sign == b.sign) {
int c = cmp(dig, b.dig);
if (c == 0) return BigInteger().purify();
else if (c == 1) return BigInteger(sub(dig, b.dig), sign).purify();
else return BigInteger(sub(b.dig, dig), (Sign)!sign).purify();
}
else {
return *this+ -b;
}
}
BigInteger operator*(const BigInteger& b) const {
BigInteger rst(DigSet(dig.size() + b.dig.size()), Sign(sign!=b.sign));
for (int i = 0; i < dig.size(); i++) {
for (int j = 0; j < b.dig.size(); j++) {
rst.dig[i + j] += dig[i] * b.dig[j];
rst.dig[i + j + 1] += rst.dig[i + j] / D;
rst.dig[i + j] %= D;
}
}
return rst.purify();
}
BigInteger operator/ (const BigInteger& b) const {
BigInteger rst(DigSet(dig.size()), Sign(sign != b.sign));
BigInteger div;
for (int i = dig.size()-1; i>=0; i--) {
div = div * D + dig[i];
BigInteger mul(b.dig, positive);
for (int j = 0x2000; j; j >>= 1) {
if (mul*j <= div) {
rst.dig[i] |= j;
div = div - mul * j;
}
}
}
return rst.purify();
}
BigInteger operator% (const BigInteger& b) const {
BigInteger rst = *this - *this / b * b;
return rst.sign == positive ? rst : b.sign == positive ? rst + b : rst - b;
}
};

int main() {
BigInteger a, b;
while (cin >> a >> b) {
cout << a + b << endl;
cout << a - b << endl;
cout << a * b << endl;
if (!(b == 0)) {
cout << a / b << endl;
cout << a % b << endl;
}
}
return 0;
}

1~n中与m互质的数的和、平方和

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
LL fact[30];
int nf;

void makeFact(LL x) {
nf = 0;
for(int i=2; (LL)i*i<=x; i++) {
if(x%i==0) {
fact[nf++] = i;
do x /= i; while (x%i == 0);
}
}
if(x>1) fact[nf++] = x;
}

inline LL sum1ToX(LL x) {
return (1+x)*x/2;
}

inline LL sqrSum1ToX(LL x) {
return x*(x+1)*(2*x+1)/6;
}

// 1~n与m互质的和(与n互质的和见欧拉函数部分)
LL sumCoprime(LL n, LL m) {
makeFact(m);
LL rst = 0;
for(int i=1; i<1<<nf; i++) {
LL f = 1;
int cnt = 0;
for(int j=0; j<nf; j++) {
if(i&(1<<j)) {
f *= fact[j];
cnt++;
}
}
if(cnt&1) rst += f*sum1ToX(n/f);
else rst -= f*sum1ToX(n/f);
}
return sum1ToX(n)-rst;
}

// 1~n与m互质的平方和
LL sqrSumCoprime(LL n, LL m) {
makeFact(m);
LL rst = 0;
for(int i=1; i<1<<nf; i++) {
LL f = 1;
int cnt = 0;
for(int j=0; j<nf; j++) {
if(i&(1<<j)) {
f *= fact[j];
cnt++;
}
}
if(cnt&1) rst += f*f*sqrSum1ToX(n/f);
else rst -= f*f*sqrSum1ToX(n/f);
}
return sqrSum1ToX(n)-rst;
}

欧拉函数

n的欧拉函数值\(\varphi(n)\)表示1~n中与n互质的数的个数。

\(n=p_1^{k_1}p_2^{k_2} \dots p_r^{k_r}\),那么: \[ \varphi(n)=\prod{p_i^{k_i-1}(p_i-1)}=n\prod(1-\frac1{p_i})=\frac{n}{\prod{p_i}}\prod{(p_i-1)} \] 另外有性质: \[ \sum_{d|n}{\varphi(d)}=n \] 经莫比乌斯反演可得: \[ \varphi(n)=\sum_{d|n}{d\cdot\mu(n/d)} \] 其他性质: \[ \sum_{d|n}{\frac{\mu^2(d)}{\varphi(d)}}=\frac{n}{\varphi(n)} \]

\[ 1到n中与n互质的数的和\sum_{\substack{1\leq k \leq n \\ (k,n)=1}}{k}=\frac{n\varphi(n)}{2} \]

1
2
3
4
5
6
7
8
9
10
11
12
LL phi(LL n) {
LL rst = n;
for (int i = 2; (LL)i*i <= n; i++) {
if (n%i == 0) {
rst -= rst / i;
do n /= i; while (n%i == 0);
}
}
if (n > 1)
rst -= rst / n;
return rst;
}

\(ax+by=c\)的整数解

\(g=gcd(a, b)\),那么当\(g\nmid c\)时,没有整数解。

先用扩展欧几里得求出\(ax'+by'=g\)的特解\(x'_0,y'_0\)

注意,原方程的通解不能表示为\(x=c/g\times x'_0, y=c/g\times y'_0\),而应该表示为

\[ \begin{cases} x=\frac cg\cdot x'_0+\frac bg\cdot t \\ y=\frac cg\cdot y'_0-\frac ag\cdot t \end{cases} \]

根据上式,如果要求\(x\)\(y\)的最小非负整数解也很容易了,即\(\frac cg x'_0\%\frac bg\)\(\frac cg y'_0\%\frac ag\)

1
2
3
4
5
6
7
8
9
10
11
12
LL xg0, yg0;
LL g = exgcd(a, b, xg0, yg0);
if(c%g) {
puts("No solution");
return 0;
}
// x为最小非负整数
LL x1 = posmod(c/g*xg0, b/g);
LL y1 = (c-a*x1)/b;
// y为最小非负整数
LL y2 = posmod(c/g*yg0, a/g);
LL x2 = (c-b*y2)/a;

Lucas定理

计算\(\binom{n}{m}\%p\),其中n和m很大而p为质数且不是很大的情况。

将n和m按p进制展开,即设: \[ n=n_kp^k+n_{k-1}p^{k-1}+\dots+n_1p+n_0 \]

\[ m=m_kp^k+m_{k-1}p^{k-1}+\dots+m_1p+m_0 \]

则: \[ \binom nm \equiv \prod \binom{n_i}{m_i} \pmod p \]

1
2
3
4
5
6
7
8
9
10
11
12
inline LL lucas(LL n, LL m, int p) {
if (n < m) return 0;
LL rst = 1;
while (m) {
int nn = n % p, mm = m % p;
if (nn < mm) return 0;
rst = rst * comb(nn, mm) % p;
n /= p;
m /= p;
}
return rst;
}

模线性方程组

解方程组: \[ x \equiv a_i \pmod{m_i} \]

\(m_i\)两两互质

\[ M=\prod{m_i} \\ W_i=\frac M{m_i} \\ x\equiv\sum{a_iW_i(W_i^{-1}\bmod m_i)}\pmod{M} \]

\(m_i\)不满足两两互质

1
2
3
4
5
6
7
8
9
10
11
12
13
LL crt(LL m[], LL r[], int n) {
LL M = m[0], R = r[0];
for (int i = 1; i < n; i++) {
LL x, y;
LL g = exgcd(M, m[i], x, y);
if ((r[i] - R) % g) return -1;
x = (r[i] - R) / g * x % (m[i] / g);
R += x * M;
M = M / g * m[i];
R %= M;
}
return (R+M)%M;
}

Splay

Splay List模板。

支持:初始化,单点增删改查,区间增删改、最值、求和、翻转、提取。

TODO: 分割,合并。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <vector>
#include <stack>
using namespace std;
typedef long long LL;
const int MAXN = 1000010;
const LL INF = 0x3f3f3f3f3f3f3f3f;

struct Node {
LL val;
int ch[2], fa, sz;
bool toSet;
LL set, add;
bool rev;
LL min, max, sum;
};
static Node V[MAXN];
static int pool[MAXN], con;

struct SplayList {
static void makePool() {
con = 0;
for (int i = 1; i < MAXN; i++) {
pool[i] = i;
}
}
static int alloc(LL val) {
int r = pool[++con];
V[r] = Node();
V[r].val = V[r].min = V[r].max = V[r].sum = val;
V[r].sz = 1;
return r;
}
static void free(int x) {
pool[con--] = x;
}
// 如果需要O(logn)的区间删除,就清空freeRec函数体
static void freeRec(int x) {
//if (V[x].ch[0]) freeRec(V[x].ch[0]);
//if (V[x].ch[1]) freeRec(V[x].ch[1]);
//free(x);
}
static void pushUp(int x) {
if (!x) return;
int& ls = V[x].ch[0], &rs = V[x].ch[1];
V[x].sz = 1;
if (ls) V[x].sz += V[ls].sz;
if (rs) V[x].sz += V[rs].sz;
V[x].min = V[x].val;
if (ls) V[x].min = min(V[x].min, V[ls].min);
if (rs) V[x].min = min(V[x].min, V[rs].min);
V[x].max = V[x].val;
if (ls) V[x].max = max(V[x].max, V[ls].max);
if (rs) V[x].max = max(V[x].max, V[rs].max);
V[x].sum = V[x].val;
if (ls) V[x].sum += V[ls].sum;
if (rs) V[x].sum += V[rs].sum;
}
// tag表示该节点已操作,但子节点未操作
static void pushDown(int x) {
if (!x) return;
int& ls = V[x].ch[0], &rs = V[x].ch[1];
if (V[x].rev) {
if (ls) {
swap(V[ls].ch[0], V[ls].ch[1]);
V[ls].rev = !V[ls].rev;
}
if (rs) {
swap(V[rs].ch[0], V[rs].ch[1]);
V[rs].rev = !V[rs].rev;
}
V[x].rev = false;
}
if (V[x].toSet) {
if (ls) {
V[ls].val = V[ls].min = V[ls].max = V[x].set;
V[ls].sum = V[x].set*V[ls].sz;
V[ls].toSet = true;
V[ls].set = V[x].set;
V[ls].add = 0;
}
if (rs) {
V[rs].val = V[rs].min = V[rs].max = V[x].set;
V[rs].sum = V[x].set*V[rs].sz;
V[rs].toSet = true;
V[rs].set = V[x].set;
V[rs].add = 0;
}
V[x].toSet = false;
}
if (V[x].add) {
if (ls) {
V[ls].val += V[x].add;
V[ls].min += V[x].add;
V[ls].max += V[x].add;
V[ls].sum += V[x].add*V[ls].sz;
V[ls].add += V[x].add;
}
if (rs) {
V[rs].val += V[x].add;
V[rs].min += V[x].add;
V[rs].max += V[x].add;
V[rs].sum += V[x].add*V[rs].sz;
V[rs].add += V[x].add;
}
V[x].add = 0;
}
}
static void rotate(int x) {
int y = V[x].fa, z = V[y].fa;
int tx = V[y].ch[1] == x, ty = V[z].ch[1] == y;
V[z].ch[ty] = x;
V[x].fa = z;
V[y].ch[tx] = V[x].ch[!tx];
if(V[x].ch[!tx]) V[V[x].ch[!tx]].fa = y;
V[x].ch[!tx] = y;
V[y].fa = x;
pushUp(y); pushUp(x);
}

int root = 0;
static int genTree(LL val[], int l, int r) {
if (l > r) return 0;
int mid = l + r >> 1;
int x = alloc(val[mid]);
V[x].ch[0] = genTree(val, l, mid - 1);
V[x].ch[1] = genTree(val, mid + 1, r);
V[V[x].ch[0]].fa = V[V[x].ch[1]].fa = x;
pushUp(x);
return x;
}
SplayList(LL val[]=nullptr, int n=0) {
root = genTree(val, 0, n - 1);
}
vector<LL> toVector() const {
vector<LL> rst(size());
stack<int> stk;
int p = root, q = 0;
while (p || !stk.empty()) {
while (p) {
pushDown(p);
stk.push(p);
p = V[p].ch[0];
}
if (!stk.empty()) {
p = stk.top();
stk.pop();
rst[q++] = V[p].val;
p = V[p].ch[1];
}
}
return rst;
}
int size() const {
return V[root].sz;
}
void splay(int x, int goal) {
pushDown(x);
while (V[x].fa != goal) {
int y = V[x].fa, z = V[y].fa;
pushDown(z); pushDown(y); pushDown(x);
int tx = V[y].ch[1] == x, ty = V[z].ch[1] == y;
if (z != goal) {
if (tx == ty) rotate(y);
else rotate(x);
}
rotate(x);
}
if (goal == 0) root = x;
}
// pos: 1~n
int at(int pos) const {
int p = root;
while (p) {
pushDown(p);
if (V[V[p].ch[0]].sz == pos - 1) return p;
else if (V[V[p].ch[0]].sz >= pos) p = V[p].ch[0];
else {
pos -= V[V[p].ch[0]].sz + 1;
p = V[p].ch[1];
}
}
return 0;
}
// pos: 0~n
int insert(int pos, LL val) {
int nn = alloc(val);
if (root == 0) {
root = nn;
}
else if (pos == 0) {
splay(at(1), 0);
V[root].ch[0] = nn;
V[nn].fa = root;
pushUp(root);
}
else if (pos == size()) {
splay(at(pos), 0);
V[root].ch[1] = nn;
V[nn].fa = root;
pushUp(root);
}
else {
splay(at(pos), 0);
int fa;
splay(fa = at(pos + 1), root);
V[fa].ch[0] = nn;
V[nn].fa = fa;
pushUp(fa);
pushUp(root);
}
return nn;
}
// pos: 1~n
void erase(int pos) {
if (pos == 1) {
splay(at(1), 0);
int x = root;
root = V[root].ch[1];
V[root].fa = 0;
free(x);
pushDown(root);
pushUp(root);
}
else if (pos == size()) {
splay(at(pos), 0);
int x = root;
root = V[root].ch[0];
V[root].fa = 0;
free(x);
pushDown(root);
pushUp(root);
}
else {
splay(at(pos - 1), 0);
int fa;
splay(fa = at(pos + 1), root);
int x = V[fa].ch[0];
V[fa].ch[0] = 0;
free(x);
pushDown(fa);
pushUp(fa);
pushDown(root);
pushUp(root);
}
}
void modify(int pos, LL val) {
splay(at(pos), 0);
V[root].val = val;
pushUp(root);
}
LL operator[](int pos) const {
return V[at(pos)].val;
}
int range(int l, int r) {
if (l == 1 && r == size()) return root;
else if (l == 1) {
splay(at(r + 1), 0);
return V[root].ch[0];
}
else if (r == size()) {
splay(at(l - 1), 0);
return V[root].ch[1];
}
else {
splay(at(l - 1), 0);
int fa;
splay(fa = at(r + 1), root);
return V[fa].ch[0];
}
}
LL Min(int l, int r) {
return V[range(l, r)].min;
}
LL Max(int l, int r) {
return V[range(l, r)].max;
}
LL sum(int l, int r) {
return V[range(l, r)].sum;
}
void reverse(int l, int r) {
int x = range(l, r);
swap(V[x].ch[0], V[x].ch[1]);
V[x].rev = !V[x].rev;
// 这里不pushup是因为其他信息不受影响
}
void modify(int l, int r, LL val) {
int x = range(l, r);
V[x].val = V[x].min = V[x].max = val;
V[x].sum = val * V[x].sz;
V[x].toSet = true;
V[x].set = val;
V[x].add = 0;
pushUp(V[x].fa);
pushUp(V[V[x].fa].fa);
}
void add(int l, int r, LL delta) {
int x = range(l, r);
V[x].val += delta;
V[x].min += delta;
V[x].max += delta;
V[x].sum += delta * V[x].sz;
V[x].add += delta;
pushUp(V[x].fa);
pushUp(V[V[x].fa].fa);
}
void erase(int l, int r) {
int x = range(l, r);
int y = V[x].fa;
if (y == 0) {
root = 0;
}
else {
int t = V[y].ch[1] == x;
V[y].ch[t] = 0;
pushUp(y);
pushUp(V[y].fa);
}
freeRec(x);
}
};

除法取整分块

相同的\(\lfloor \frac ni \rfloor\)分为一块。

1
2
3
4
5
6
7
8
9
10
#include <cstdio>

int main() {
int n = 50;
for (int i = 1, j=0; i <= n; i = j + 1) {
int temp = j;
j = n / (n / i);
printf("cnt:%d %d/(%d..%d)=%d\n", j-temp, n, temp+1, j, n/j);
}
}

Heron法开方(其实也是Newton法)

返回平方根向下取整,较快。

1
2
3
4
5
6
7
static BigInteger sqrt(BigInteger n) {
BigInteger rst = BigInteger.ONE.shiftLeft(n.bitLength()>>1);
while(rst.pow(2).compareTo(n)>0 || rst.add(BigInteger.ONE).pow(2).compareTo(n)<=0) {
rst = n.divide(rst).add(rst).shiftRight(1);
}
return rst;
}

Newton法开方

整数

较慢

1
2
3
4
5
6
7
8
static BigInteger sqrt(BigInteger n) {
BigInteger x = n, y = n.add(BigInteger.ONE).shiftRight(1);
while(y.compareTo(x)<0) {
x = y;
y = x.add(n.divide(x)).shiftRight(1);
}
return x;
}

小数

推荐。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
static final BigDecimal TWO = BigDecimal.valueOf(2);
static BigDecimal sqrt(BigDecimal n, final int SCALE) {
double dn = n.doubleValue();
BigDecimal temp = null, x;
if(Double.isFinite(dn)) {
if(Math.abs(dn)<1e-10) return ZERO;
x = BigDecimal.valueOf(Math.sqrt(n.doubleValue()));
}
else
x = new BigDecimal(BigInteger.ONE.shiftLeft(n.toBigInteger().bitLength()>>1));
while(!x.equals(temp)) {
temp = x;
x = n.divide(x, SCALE, ROUND_HALF_EVEN).add(x).divide(TWO, SCALE, ROUND_HALF_EVEN);
}
return x;
}

平方剩余

模数为奇质数

判断\(x^2\equiv a \pmod p\)是否有解,要求p是奇质数。

用多次测试的方法可以检测一个大数是否为完全平方数。

1
2
3
4
inline bool quadraticResidue(LL a, LL p) {
a%=p;
return a==0 || fastPow(a, p >> 1, p) == 1;
}

模数为奇数

分解\(m=p_1^{k_1}p_2^{k_2} \dots p_r^{k_r}\),然后逐个判断quadraticResidue(a, p_i)

线段树

区间加、乘。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
/// 洛谷 P3373
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std;
typedef long long LL;

const int MAXN = 1e5 + 10;
int n, m;
LL mod;

struct {
LL sum;
LL add, mul;
int l, r;
int len() const {
return r - l + 1;
}
} T[MAXN*4];
LL a[MAXN];

#define LS (O<<1)
#define RS (O<<1|1)
#define M (T[O].l+T[O].r>>1)

inline void pullUp(int O) {
T[O].sum = (T[LS].sum + T[RS].sum) % mod;
}

inline void pushDown(int O) {
LL& add = T[O].add, &mul = T[O].mul;
T[LS].sum = T[LS].sum*mul%mod;
T[LS].mul = T[LS].mul*mul%mod;
T[LS].add = T[LS].add*mul%mod;
T[RS].sum = T[RS].sum*mul%mod;
T[RS].mul = T[RS].mul*mul%mod;
T[RS].add = T[RS].add*mul%mod;
mul = 1;
T[LS].sum = (T[LS].sum + add*T[LS].len()%mod) % mod;
T[LS].add = (T[LS].add + add) % mod;
T[RS].sum = (T[RS].sum + add*T[RS].len()%mod) % mod;
T[RS].add = (T[RS].add + add) % mod;
add = 0;
}

inline LL query(int l, int r, int O) {
if (l <= T[O].l && T[O].r <= r) {
return T[O].sum;
}
pushDown(O);
LL rst = 0;
if (l <= M)
rst = query(l, r, LS);
if (r > M)
rst = (rst + query(l, r, RS)) % mod;
pullUp(O);
return rst;
}

inline void updateAdd(int l, int r, LL d, int O) {
if (l <= T[O].l && T[O].r <= r) {
T[O].sum = (T[O].sum + d * T[O].len() % mod) % mod;
T[O].add = (T[O].add + d) % mod;
return;
}
pushDown(O);
if (l <= M)
updateAdd(l, r, d, LS);
if (r > M)
updateAdd(l, r, d, RS);
pullUp(O);
}

inline void updateMul(int l, int r, LL d, int O) {
if (l <= T[O].l && T[O].r <= r) {
T[O].sum = T[O].sum*d % mod;
T[O].add = T[O].add*d % mod;
T[O].mul = T[O].mul*d % mod;
return;
}
pushDown(O);
if (l <= M)
updateMul(l, r, d, LS);
if (r > M)
updateMul(l, r, d, RS);
pullUp(O);
}

inline void build(int O, int L, int R) {
T[O] = { 0,0,1,L,R };
if (L == R) {
T[O].sum = a[L];
return;
}
build(LS, L, M);
build(RS, M+1, R);
pullUp(O);
}

int main() {
while (scanf("%d%d%lld", &n, &m, &mod) == 3) {
for (int i = 1; i <= n; i++) {
scanf("%lld", a + i);
}
build(1, 1, n);
for (int i = 0; i < m; i++) {
int op;
scanf("%d", &op);
if (op == 1) {
int x, y;
LL k;
scanf("%d%d%lld", &x, &y, &k);
updateMul(x, y, k, 1);
}
else if (op == 2) {
int x, y;
LL k;
scanf("%d%d%lld", &x, &y, &k);
updateAdd(x, y, k, 1);
}
else if (op == 3) {
int x, y;
scanf("%d%d", &x, &y);
printf("%lld\n", query(x, y, 1));
}
}
}
return 0;
}

后缀自动机

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
/// hihoCoder 1445
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std;

const int MAXL = 1e6+10, MAXHA = 26;
struct SuffixAutomaton {
struct {
int par, ch[MAXHA];
int maxlen;
} v[MAXL*2];

int cos = 1, s = 1, last = 1;
int newState() {
return ++cos;
}
int push(int ha) {
int p = last;
const int np = newState();
last = np;
v[np].maxlen = v[p].maxlen + 1;
while (p && !v[p].ch[ha]) {
v[p].ch[ha] = np;
p = v[p].par;
}
if (!p) {
v[np].par = s;
return np;
}
const int q = v[p].ch[ha];
if (v[q].maxlen == v[p].maxlen + 1)
v[np].par = q;
else {
const int nq = newState();
v[nq].maxlen = v[p].maxlen + 1;
memcpy(v[nq].ch, v[q].ch, sizeof v[q].ch);
v[nq].par = v[q].par;
v[q].par = v[np].par = nq;
while (v[p].ch[ha] == q) {
v[p].ch[ha] = nq;
p = v[p].par;
}
}
return np;
}
};

char a[MAXL];
SuffixAutomaton sa;

int main() {
scanf("%s", a);
for (char* p = a; *p; ++p) {
sa.push(*p - 'a');
}
long long ans = 0;
for (int i = 2; i <= sa.cos; i++) {
ans += sa.v[i].maxlen - sa.v[sa.v[i].par].maxlen;
}
printf("%lld\n", ans);
return 0;
}

Tarjan强连通分量

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
/// HDU 1269
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std;

const int MAXN = 1e5 + 10, MAXM = 1e6 + 10;
int n, m;
struct {
int v, next;
} e[MAXM];
int head[MAXN], coe;

void addEdge(int u, int v) {
e[coe].v = v;
e[coe].next = head[u];
head[u] = coe++;
}

int idx, stk[MAXN], ps, dfn[MAXN], low[MAXN];
bool instack[MAXN];
int belong[MAXN], nscc;

void tarjan(int u) {
dfn[u] = low[u] = ++idx;
stk[ps++] = u;
instack[u] = true;
for (int i = head[u]; i != -1; i = e[i].next) {
int v = e[i].v;
if (!dfn[v]) {
tarjan(v);
low[u] = min(low[u], low[v]);
}
else if (instack[v]) {
low[u] = min(low[u], dfn[v]);
}
}
if (low[u] == dfn[u]) {
++nscc;
int temp;
do {
temp = stk[--ps];
belong[temp] = nscc;
instack[temp] = false;
} while (temp != u);
}
}

void scc() {
idx = 0;
nscc = 0;
memset(dfn, 0, sizeof dfn);
for (int i = 1; i <= n; i++) {
if (!dfn[i])
tarjan(i);
}
}

int main() {
while (scanf("%d%d", &n, &m) == 2 && (n||m)) {
memset(head, -1, sizeof head);
coe = 0;
for (int i = 0; i < m; i++) {
int a, b;
scanf("%d%d", &a, &b);
addEdge(a, b);
}
scc();
puts(nscc == 1 ? "Yes" : "No");
}
return 0;
}

Tarjan割点与桥

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
int idx, low[MAXN], dfn[MAXN];
int ncc[MAXN]; // 删除当前点,所在联通分量会变为几块。大于1是割点,等于0是孤立点

// 对每个连通块执行一次tarjan(u, -1)
void tarjan(int u, int fa) {
dfn[u] = low[u] = ++idx;
ncc[u] = fa != -1;
for (int i = head[u]; i != -1; i = e[i].next) {
int v = e[i].v;
if (v == fa) continue;
if (!dfn[v]) {
tarjan(v, u);
low[u] = min(low[u], low[v]);
if (fa==-1 || dfn[u] <= low[v])
ncc[u]++;
if (dfn[u] < low[v]) {
// uv是桥
}
}
else {
low[u] = min(low[u], dfn[v]);
}
}
}

树分治

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
/// POJ 1741
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <utility>
#include <vector>
using namespace std;
typedef long long LL;
const int MAXN = 10010;
int n;
LL k;
struct {
int v, next;
LL w;
} e[MAXN*2];
int head[MAXN], coe;

void addEdge(int u, int v, LL w) {
e[coe].v = v;
e[coe].w = w;
e[coe].next = head[u];
head[u] = coe++;
}

bool rooted[MAXN];
int sz[MAXN], maxChildSz[MAXN];

int dfsSize(int u, int fa) {
sz[u] = 1;
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
if(v==fa || rooted[v]) continue;
sz[u] += dfsSize(v, u);
}
return sz[u];
}

void dfsMaxChild(int u, int fa, int r) {
maxChildSz[u] = sz[r]-sz[u];
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
if(v==fa || rooted[v]) continue;
maxChildSz[u] = max(maxChildSz[u], sz[v]);
dfsMaxChild(v, u, r);
}
}

int dfsRoot(int u, int fa) {
int rst = u;
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
if(v==fa || rooted[v]) continue;
int t = dfsRoot(v, u);
if(maxChildSz[t]<maxChildSz[rst]) {
rst = t;
}
}
return rst;
}

vector<int> dis;

void dfsDis(int u, int fa, int d) {
dis.push_back(d);
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
if(v==fa || rooted[v]) continue;
dfsDis(v, u, d+e[i].w);
}
}

int cal(int u, int d) {
int rst = 0;
dis.clear();
dfsDis(u, -1, d);
sort(dis.begin(), dis.end());
int i=0, j=dis.size()-1;
while(i<j) {
while(dis[i]+dis[j]>k && i<j) j--;
rst += j-i;
i++;
}
return rst;
}

int ans;
void treeDivide(int u) {
dfsSize(u, -1);
dfsMaxChild(u, -1, u);
u = dfsRoot(u, -1);
ans += cal(u, 0);
rooted[u] = true;
for(int i=head[u]; i!=-1; i=e[i].next) {
int v = e[i].v;
if(rooted[v]) continue;
ans -= cal(v, e[i].w);
treeDivide(v);
}
}

int main() {
while(scanf("%d%lld", &n, &k)==2 && (n||k)) {
memset(head, -1, sizeof head);
memset(rooted, 0, sizeof rooted);
coe = 0;
for(int i=1; i<n; i++) {
int a,b,c;
scanf("%d%d%d",&a, &b, &c);
addEdge(a,b,c);
addEdge(b,a,c);
}
ans = 0;
treeDivide(1);
printf("%d\n", ans);
}
return 0;
}

KD树

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
/// HDU 4347
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <queue>
using namespace std;
typedef long long LL;
const int MAXN = 200010, MAXK = 5;

int n,k;

struct Point {
int x[MAXK];
int& operator[] (int i) {
return x[i];
}
};
LL dis2(Point& a, Point& b) {
LL rst = 0;
for(int i=0; i<k; i++) {
rst += (LL)(a[i]-b[i])*(a[i]-b[i]);
}
return rst;
}

struct Node {
Point x;
int ch[2];
} node[MAXN];

struct KDTree {
int root;

inline int build(int l, int r, int d) {
if(l>r) return 0;
int mid = l+r>>1;
nth_element(node+l, node+mid, node+r+1, [d](Node& a, Node& b) {
return a.x[d%k]<b.x[d%k];
});
node[mid].ch[0] = build(l, mid-1, d+1);
node[mid].ch[1] = build(mid+1, r, d+1);
return mid;
}
void build(int len) {
root = build(1, len, 0);
}

inline void query(int p, Point x, priority_queue<pair<LL,int> >& que, int nk, int d) {
pair<LL,int> now = make_pair(dis2(x, node[p].x), p);
que.push(now);
if(que.size()>nk)
que.pop();
LL t = node[p].x[d%k]-x[d%k];
int a = node[p].ch[0], b = node[p].ch[1];
if(x[d%k]>=node[p].x[d%k])
swap(a, b);
if(a)
query(a, x, que, nk, d+1);
if(b && que.top().first>t*t)
query(b, x, que, nk, d+1);
}
vector<Point> knn(Point x, int nk) {
priority_queue<pair<LL,int> > que; // dis,id
query(root, x, que, nk, 0);
vector<Point> rst(que.size());
for(int i=que.size()-1; i>=0; i--) {
rst[i] = node[que.top().second].x;
que.pop();
}
return rst;
}
} tree;

int main() {
while(scanf("%d%d", &n, &k)==2) {
for(int i=1; i<=n; i++) {
for(int j=0; j<k; j++) {
scanf("%d", &node[i].x[j]);
}
}
tree.build(n);
int q;
scanf("%d", &q);
for(int i=0; i<q; i++) {
Point x;
int m;
for(int j=0; j<k; j++) {
scanf("%d", &x[j]);
}
scanf("%d", &m);
vector<Point> rst = tree.knn(x, m);
printf("the closest %d points are:\n", m);
for(Point& it : rst) {
for(int j=0; j<k; j++) {
printf("%d%c", it[j], " \n"[j==k-1]);
}
}
}
}
return 0;
}

树链剖分

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
int id[MAXN],   //树上结点在线段树的编号
aid[MAXN], //线段树结点在树上的编号
son[MAXN], //重儿子
dep[MAXN], fa[MAXN],
siz[MAXN], //子树大小
top[MAXN], //所在重链的根节点
tot; // dfs计时器

// 先调用dfs1(1, 0)
void dfs1(int x, int f) {
fa[x] = f;
dep[x] = dep[f] + 1;
son[x] = 0;
siz[x] = 1;
for (int i = head[x]; i != -1; i = e[i].next) {
int v = e[i].v;
if (v == f) continue;
dfs1(v, x);
siz[x] += siz[v];
if (!son[x] || siz[son[x]] < siz[v]) son[x] = v;
}
}

// 再调用dfs2(1, 1)
void dfs2(int x, int tp) {
top[x] = tp;
id[x] = ++tot;
aid[tot] = x;
if (son[x]) dfs2(son[x], tp);
for (int i = head[x]; i != -1; i = e[i].next) {
int v = e[i].v;
if (v != fa[x] && v != son[x]) {
dfs2(v, v);
}
}
}

// 查询u到v
int qsum(int u, int v) {
int rst = 0;
while (top[u] != top[v]) {
if (dep[top[u]] < dep[top[v]]) std::swap(u, v);
// 处理线性区间 id[top[u]] -> id[u]
rst += querySum(1, 1, n, id[top[u]], id[u]);
u = fa[top[u]];
}
if (dep[u] > dep[v]) std::swap(u, v);
// 最后处理线性区间 id[u] -> id[v]
return rst + querySum(1, 1, n, id[u], id[v]);
}

树状数组

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
LL C[MAXN], CI[MAXN];

// 给l~n每个数字加d
void add(int l, LL d) {
for (int i = l; i <= n; i+=i&-i) {
C[i] += d;
CI[i] += l * d;
}
}

void add(int l, int r, LL d) {
add(l, d);
add(r + 1, -d);
}

// 查询1~r的和
LL sum(int r) {
LL rst = 0;
for (int i = r; i; i ^= i & -i) {
rst += (r + 1)*C[i] - CI[i];
}
return rst;
}

最大密度子图

选择一个子图,使得边数/点数最大。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
/// UVALive 7037
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
#include <queue>
#include <utility>
using namespace std;
// double型的最大流
const int MAXN = 1e6 + 10, MAXM = 1e6 + 10;
const double INF = 1e50;
const double EPS = 1e-6;

struct Edge {
int v, next;
double cap, flow, cost;
} e[MAXM];
int head[MAXN], coe;
int N;

inline int dcmp(double x, double y) {
if(x>y+EPS) return 1;
if(x+EPS<y) return -1;
return 0;
}

// 每次初始化,节点编号为0~n-1
void init(int n) {
N = n;
coe = 0;
memset(head, -1, sizeof(int)*n);
}
void addEdge(int u, int v, double cap, double cost) {
e[coe].v = v;
e[coe].cap = cap;
e[coe].cost = cost;
e[coe].flow = 0;
e[coe].next = head[u];
head[u] = coe++;
e[coe].v = u;
e[coe].cap = 0;
e[coe].cost = -cost;
e[coe].flow = 0;
e[coe].next = head[v];
head[v] = coe++;
}

int dep[MAXN];

bool dinic_bfs(int source, int dest) {
memset(dep, -1, sizeof(int)*N);
dep[source] = 0;
queue<int> que;
que.push(source);
while (!que.empty()) {
int i = que.front();
que.pop();
for (int j = head[i]; j != -1; j = e[j].next) {
if (dep[e[j].v] < 0 && dcmp(e[j].cap,e[j].flow)>0) {
dep[e[j].v] = dep[i] + 1;
que.push(e[j].v);
}
}
}
return dep[dest] > 0;
}

inline double dinic_find(int x, double low, int source, int dest) {
if (dcmp(low, 0)<=0) return false;
if (x == dest) return low;
double cost = 0;
for (int i = head[x]; i != -1; i = e[i].next) {
if ( dcmp(e[i].cap, e[i].flow)>0 && dep[e[i].v] == dep[x] + 1) {
double a = dinic_find(e[i].v, min(low - cost, e[i].cap-e[i].flow), source, dest);
if (dcmp(a, 0)>0) {
cost += a;
e[i].flow += a;
e[i ^ 1].flow -= a;
if (dcmp(cost, low)>=0)
break;
}
else {
dep[e[i].v] = -1;
}
}
}
return cost;
}

double dinic(int source, int dest) {
double ans = 0;
while (dinic_bfs(source, dest)) {
double tans;
while (dcmp(tans = dinic_find(source, INF, source, dest), 0)>0)
ans += tans;
}
return ans;
}

int n;
int a[MAXN];
int deg[MAXN];
Edge backupE[MAXM];
int backupHead[MAXN], backupCoe;

int main() {
int T;
scanf("%d", &T);
for(int cs=1; cs<=T; cs++) {
scanf("%d", &n);
memset(deg, 0, sizeof(int)*(n+2));
for(int i=1; i<=n; i++) {
scanf("%d", a+i);
}
int source = 0, dest = n+1;
init(dest+1);
int m = 0;
for(int i=1; i<=n; i++) {
for(int j=i+1; j<=n; j++) {
if(a[i]>a[j]) {
// 一定是无向图,记录边数和度数
m++;
deg[i]++;
deg[j]++;
addEdge(i, j, 1, 0);
addEdge(j, i, 1, 0);
}
}
}
for(int i=1; i<=n; i++) {
addEdge(source, i, m, 0); // 首先连接source->每个点,容量为边数
}

memcpy(backupE, e, sizeof(Edge)*coe);
memcpy(backupHead, head, sizeof(int)*N);
backupCoe = coe;

double left = 0, right = m;
while(left+1e-9<right) {
double mid = (left+right)/2;
memcpy(e, backupE, sizeof(Edge)*coe);
memcpy(head, backupHead, sizeof(int)*N);
coe = backupCoe;
for(int i=1; i<=n; i++) {
addEdge(i, dest, m+2*mid-deg[i], 0); // 每次连接每个点->dest,容量为 边数+2*mid-度数
}
double maxflow = dinic(source, dest);
// 二分策略
if(n*m>maxflow) {
left = mid;
}
else {
right = mid;
}
}
printf("Case #%d: %.8f\n", cs, (left+right)/2);
}
return 0;
}

最大权闭合子图

有一个有向图,每一个点都有一个权值(可以为正或负或0),选择一个权值和最大的子图,使得每个点的后继都在子图里面,这个子图就叫最大权闭合子图。

从源点s向每个正权点连一条容量为权值的边,每个负权点向汇点t连一条容量为权值的绝对值的边,有向图原来的边容量全部为无限大。

求它的最小割,割掉后,与源点s连通的点构成最大权闭合子图,权值为正权值之和-最小割。

约瑟夫环

编号为0~n-1

1
2
f(1)=0
f(i)=(f(i-1)+k)%i

优化做法

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
LL josephus(LL n, LL k) {
if (k == 1) return n - 1;
LL ans = 0;
for (LL i = 2; i <= n;) {
if (ans + k >= i) {
ans = (ans + k) % i;
i++;
continue;
}
LL step = (i - ans - 2) / (k - 1);
if (i + step > n) {
ans += (n - (i - 1))*k;
break;
}
i += step; ans += step * k;
}
return ans%n;
}

博弈论

Bash游戏

只有一堆n个物品,两个人轮流从这堆物品中取物,规定每次至少取1个,最多取m个。最后取光者得胜。

只要n%(m+1)!=0,则先取者一定获胜。

Nim游戏

有若干堆石子,每堆石子的数量都是有限的,合法的移动是“选择一堆石子并拿走若干颗(不能不拿)”,如果轮到某个人时所有的石子堆都已经被拿空了,则判负(因为他此刻没有任何合法的移动)。

对于一个Nim游戏的局面,它是P-position(先手必败)当且仅当a1^a2^...^an=0

Wythoff游戏

两堆石子,博弈双方每次可以取一堆石子中的任意个,不能不取,或者取两堆石子中的相同个。先取完者赢。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
import java.math.BigDecimal;
import java.math.RoundingMode;
import java.util.Scanner;

public class Main {
static final BigDecimal SQRT5 = new BigDecimal(
"2.2360679774997896964091736687312762354406183596115257242708972454105209256378048994144144083787822750");
static Scanner in = new Scanner(System.in);

public static void main(String[] args) {
while (in.hasNextBigDecimal()) {
BigDecimal a = in.nextBigDecimal(), b = in.nextBigDecimal();
BigDecimal diff = a.subtract(b).abs();
if (a.min(b).equals(diff.multiply(SQRT5.add(BigDecimal.ONE)).divide(BigDecimal.valueOf(2), 0, RoundingMode.FLOOR))) {
// min(a,b) == floor(|a-b|*(sqrt(5)+1)/2)
// 先手败
System.out.println(0);
} else {
System.out.println(1);
}
}
}
}

01分数规划

求:

\[ max\{\frac{\sum_{i\ is\ chosen}{a_i}}{\sum_{i\ is\ chosen}{b_i}}\} \] 发现: \[ \frac{\sum{a_i}}{\sum{b_i}} \geq x \Leftrightarrow \sum{a_i}-x \sum{b_i} \geq 0 \] 所以二分\(x\),对\(a_i-xb_i\)排序选前几个就可以了。

表达式处理

中缀表达式直接求值

  • 如果c='+'或者'-',那么当符号栈非空并且栈顶元素不为'('时,计算一次,计算完成后再把符号入栈。
  • 如果c='*'或者'/',那么当符号栈非空并且栈顶元素为'*'或者'/'时,计算一次,计算完成后把符号入栈。
  • 如果c是数字,那么入数字栈。
  • 如果c是'(',那么直接入符号栈。
  • 如果c是')',那么一直计算到栈顶为'('为止,最后再把'('弹出栈。

不过需要注意的是在求值之前需要对表达式进行预处理,去掉空格、识别负号(区分'-'是作为减号还是负号),提取操作数等。

对于'-'的区分,主要判别方法为:

  • 若前一个字符为'(',则必定为负号。
  • 若前一个字符为')'或者数字,则必定为减号。
  • 若前面一个字符为其他运算符,如*,/,则必定是负号。
  • 若前面没有字符,即该字符为表达式的第一个字符,则必定是负号。

也就是说只有一种情况下,'-'是作为减号使用的,就是前一个字符为')'或者数字的时候。

中缀表达式转后缀表达式

对输入的中缀表达式从左到右遍历:

  • 如果遇到数字,直接添加到后缀表达式末尾。

  • 如果遇到运算符+、-、*、/:

    先判断栈是否为空。若是,则直接将此运算符压入栈。若不是,则查看当前栈顶元素。若栈顶元素优先级大于或等于此操作符级别,则弹出栈顶元素,将栈顶元素添加到后缀表达式中,并继续进行上述判断。如果不满足上述判断或者栈为空,将这个运算符入栈。要注意的是,经过上述步骤,这个运算符最终一定会入栈。

  • 如果遇到括号:

    如果是左括号,直接入栈。如果是右括号,弹出栈中第一个左括号前所有的操作符,并将左括号弹出。(右括号别入栈)。

  • 字符串遍历结束后,如果栈不为空,则弹出栈中所有元素,将它们添加到后缀表达式的末尾,直到栈为空。

Miller-Rabin素性测试

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
const LL test[12] = { 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31,37 };

LL fastPow(__int128 base, LL exp, LL mod) {
LL rst = 1;
while (exp) {
if (exp & 1) rst = rst * base%mod;
base = base * base%mod;
exp >>= 1;
}
return rst;
}

bool millerrabin(LL n) {
LL d = n - 1;
int s = 0;
while (!(d & 1)) {
d >>= 1;
s++;
}
for (int i = 0; i < 12 && test[i]<n; i++) {
LL a = test[i];
LL t = fastPow(a, d, n);
if (t == 1)
continue;
int j;
for (j = 0; j < s; j++) {
if (t == n - 1)
break;
t = (__int128)t*t%n;
}
if (j == s)
return false;
}
return true;
}

std::regex

std::regex_match匹配整个字符串。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
#include <iostream>
#include <string>
#include <regex>

int main() {
// 简单正则表达式匹配
std::string fnames[] = { "foo.txt", "bar.txt", "baz.dat", "zoidberg" };
std::regex txt_regex("[a-z]+\\.txt");

for (const auto &fname : fnames) {
std::cout << fname << ": " << std::regex_match(fname, txt_regex) << '\n';
}
/*
foo.txt: 1
bar.txt: 1
baz.dat: 0
zoidberg: 0
*/

// 提取子匹配
std::regex base_regex("([a-z]+)\\.txt");
std::smatch base_match;

for (const auto &fname : fnames) {
if (std::regex_match(fname, base_match, base_regex)) {
// 首个 sub_match 是整个字符串;下个
// sub_match 是首个有括号表达式。
if (base_match.size() == 2) {
std::ssub_match base_sub_match = base_match[1];
std::string base = base_sub_match.str();
std::cout << fname << " has a base of " << base << '\n';
}
}
}
/*
foo.txt has a base of foo
bar.txt has a base of bar
*/

// 提取几个子匹配
std::regex pieces_regex("([a-z]+)\\.([a-z]+)");
std::smatch pieces_match;

for (const auto &fname : fnames) {
if (std::regex_match(fname, pieces_match, pieces_regex)) {
std::cout << fname << '\n';
for (size_t i = 0; i < pieces_match.size(); ++i) {
std::ssub_match sub_match = pieces_match[i];
std::string piece = sub_match.str();
std::cout << " submatch " << i << ": " << piece << '\n';
}
}
}
/*
foo.txt
submatch 0: foo.txt
submatch 1: foo
submatch 2: txt
bar.txt
submatch 0: bar.txt
submatch 1: bar
submatch 2: txt
baz.dat
submatch 0: baz.dat
submatch 1: baz
submatch 2: dat
*/
}

std::regex_search用于搜索第一个匹配的部分。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
#include <iostream>
#include <string>
#include <regex>

int main() {
std::string lines[] = { "Roses are #ff0000",
"violets are #0000ff",
"all of my base are belong to you" };

std::regex color_regex("#([a-f0-9]{2})"
"([a-f0-9]{2})"
"([a-f0-9]{2})");

// 简单匹配
for (const auto &line : lines) {
std::cout << line << ": " << std::boolalpha
<< std::regex_search(line, color_regex) << '\n';
}
std::cout << '\n';
/*
Roses are #ff0000: true
violets are #0000ff: true
all of my base are belong to you: false
*/

// 展示每个匹配中有标记子表达式的内容
std::smatch color_match;
for (const auto& line : lines) {
if (std::regex_search(line, color_match, color_regex)) {
std::cout << "matches for '" << line << "'\n";
std::cout << "Prefix: '" << color_match.prefix() << "'\n";
for (size_t i = 0; i < color_match.size(); ++i)
std::cout << i << ": " << color_match[i] << '\n';
std::cout << "Suffix: '" << color_match.suffix() << "\'\n\n";
}
}
/*
matches for 'Roses are #ff0000'
Prefix: 'Roses are '
0: #ff0000
1: ff
2: 00
3: 00
Suffix: ''

matches for 'violets are #0000ff'
Prefix: 'violets are '
0: #0000ff
1: 00
2: 00
3: ff
Suffix: ''
*/

// 重复搜索(参阅 std::regex_iterator )
std::string log(R"(
Speed: 366
Mass: 35
Speed: 378
Mass: 32
Speed: 400
Mass: 30)");
std::regex r(R"(Speed:\t\d*)");
std::smatch sm;
while (regex_search(log, sm, r)) {
std::cout << sm.str() << '\n';
log = sm.suffix();
}
/*
Speed: 366
Speed: 378
Speed: 400
*/

// C 风格字符串演示
std::cmatch cm;
if (std::regex_search("this is a test", cm, std::regex("test")))
std::cout << "\nFound " << cm[0] << " at position " << cm.prefix().length();
/*
Found test at position 10
*/
}

std::regex_iteratorstd::regex_token_iterator是匹配迭代器。重要。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
#include <cstdio>
#include <cstring>
#include <cmath>
#include <string>
#include <regex>
#include <iterator>
#include <vector>
using namespace std;

vector<string> split(const string& str, const regex& regex) {
vector<string> rst;
sregex_token_iterator it(str.begin(), str.end(), regex, -1);
sregex_token_iterator end;
while (it != end) {
rst.push_back(*it++);
}
return rst;
}

vector<string> extract(const string& str, const regex& regex) {
vector<string> rst;
sregex_iterator it(str.begin(), str.end(), regex);
sregex_iterator end;
for (; it != end; ++it) {
rst.push_back(it->str());
}
return rst;
}

int main() {
string hosts = "example.com, hello.world.com, blog.yuki-nagato.com, some.long.hostname.com";
printf("%s\n", hosts.c_str());

vector<string> hosts1 = split(hosts, regex(",\\s*"));
for (const string& host : hosts1) {
printf("host: \"%s\"\n", host.c_str());
}
putchar('\n');

vector<string> hosts2 = extract(hosts, regex("[\\w\\.-]+"));
for (const string& host : hosts2) {
printf("host: \"%s\"\n", host.c_str());
}
return 0;
}

std::regex_replace比较简单。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
#include <string>
#include <regex>
#include <iostream>
using namespace std;

int main() {
string subject("its all about geeksforgeeks");

string result1, result2, result3, result4;
string result5;

// regex object
regex re("(geeks)(.*)");

// $2 contains, 2nd capturing group which is (.*) means
// string after "geeks" which is "forgeeks". hence
// the match(geeksforgeeks) will be replaced by "forgeeks".
// so the result1 = "its all about forgeeks"
result1 = regex_replace(subject, re, "$2");

// similarly $1 contains, 1 st capturing group which is
// "geeks" so the match(geeksforgeeks) will be replaced
// by "geeks".so the result2 = "its all about geeks".
result2 = regex_replace(subject, re, "$1");

// $0 contains the whole match
// so result3 will remain same.
result3 = regex_replace(subject, re, "$0");

// $0 and $& contains the whole match
// so result3 will remain same
result4 = regex_replace(subject, re, "$&");

// Here number of capturing group
// is 2 so anything above 2
// will be replaced by nothing.
result5 = regex_replace(subject, re, "$6");

cout << result1 << endl << result2 << endl;
cout << result3 << endl << result4 << endl
<< result5;
/*
its all about forgeeks
its all about geeks
its all about $0
its all about geeksforgeeks
its all about
*/
return 0;
}
分享到 评论