算法竞赛模板
算法竞赛常用代码模板集合，包括数据结构、图论、数学和字符串算法实现。
算法竞赛模板
DataStructure
FenwickTree(Range Add)
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
template<class T>
struct Fenwick {
int n;
std::vector<T> c1, c2;
Fenwick(int n = 0) {
init(n + 1);
}
void init(int n) {
this->n = n;
c1.assign(n, T());
c2.assign(n, T());
}
void add(int p, T x) {
for (int i = p; i <= n; i += i & -i) {
c1[i] += x;
c2[i] += x * p;
}
}
void add(int l, int r, T x) {
add(l, x);
add(r + 1, -x);
}
T ask(int p) {
auto ans = T();
for (int i = p; i > 0; i -= i & -i) {
ans += (p + 1) * c1[i] - c2[i];
}
return ans;
}
T ask(int l, int r) {
return ask(r) - ask(l - 1);
}
};
FenwickTree(Single Add)
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
template <class T>
struct fenwick_tree {
public:
fenwick_tree() : _n(0) {}
explicit fenwick_tree(int n) : _n(n), data(n, 0) {}
void add(int p, T x) {
assert(0 <= p && p < _n);
p++;
while (p <= _n) {
data[p - 1] += x;
p += p & -p;
}
}
// sum of [l, r)
T sum(int l, int r) {
assert(0 <= l && l <= r && r <= _n);
return sum(r) - sum(l);
}
private:
int _n;
vector<T> data;
T sum(int r) {
T s = 0;
while (r > 0) {
s += data[r - 1];
r -= r & -r;
}
return s;
}
};
St Table
1
2
3
4
5
6
7
8
9
for (int j = 1; j <= __lg(n); j++) {
for (int i = 1; i + (1 << j) - 1 <= n; i++) {
st[i][j] = max(st[i][j - 1], st[i + (1 << j - 1)][j - 1]);
}
}
auto get = [&](int l, int r) {
int k = __lg(r - l + 1);
return max(st[l][k], st[r - (1 << k) + 1][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
struct Seg {
int l, r;
int cnt;
}tr[N << 5];
int build(int l, int r) {
int p = ++ idx;
if (l == r) {
return p;
}
int mid = l + r >> 1;
tr[p].l = build(l, mid), tr[p].r = build(mid + 1, r);
return p;
}
int insert(int p, int l, int r, int x) {
int q = ++ idx;
tr[q] = tr[p];
if (l == r) {
tr[q].cnt ++;
return q;
}
int mid = l + r >> 1;
if (x <= mid) {
tr[q].l = insert(tr[p].l, l, mid, x);
} else {
tr[q].r = insert(tr[p].r, mid + 1, r, x);
}
pushup(u);
return q;
}
int query(int q, int p, int l, int r, int k) {
if (l == r) {
return r;
}
int cnt = tr[tr[q].l].cnt - tr[tr[p].l].cnt;
int mid = l + r >> 1;
if (k <= cnt) {
return query(tr[q].l, tr[p].l, l, mid, k);
} else {
return query(tr[q].r, tr[p].r, mid + 1, r, k - cnt);
}
}
Segment Tree 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
template<class Info, class Tag>
struct LazySegmentTree {
int n;
vector<Info> info;
vector<Tag> tag;
LazySegmentTree(int n, Info v = Info()) {
init(vector(n, v));
}
template<class T>
LazySegmentTree(vector<T> a) {
init(a);
}
void pushup(int u) {
info[u] = info[u << 1] + info[u << 1 | 1];
}
void change(int u, const Tag &v) {
info[u].apply(v);
tag[u].apply(v);
}
void pushdown(int u) {
change(u << 1, tag[u]);
change(u << 1 | 1, tag[u]);
tag[u] = Tag();
}
template<class T>
void init(vector<T> a) {
n = a.size();
info.assign(n << 2, Info());
tag.assign(n << 2, Tag());
n--;
function<void(int, int, int)> build = [&](int u, int l, int r) {
if (l == r) {
info[u] = a[l];
return;
}
int mid = l + r >> 1;
build(u << 1, l, mid);
build(u << 1 | 1, mid + 1, r);
pushup(u);
};
build(1, 0, n);
}
void modify(int u, int l, int r, int x, const Info &v) {
if (l == r) {
info[u] = v;
return;
}
pushdown(u);
int mid = l + r >> 1;
if (mid >= x) {
modify(u << 1, l, mid, x, v);
} else {
modify(u << 1 | 1, mid + 1, r, x, v);
}
pushup(u);
}
void modify(int x, const Info &v) {
modify(1, 0, n, x, v);
}
void rangeModify(int u, int l, int r, int L, int R, const Tag &v) {
if (l >= L && r <= R) {
change(u, v);
return;
}
int mid = l + r >> 1;
pushdown(u);
if (mid >= L) {
rangeModify(u << 1, l, mid, L, R, v);
}
if (mid < R) {
rangeModify(u << 1 | 1, mid + 1, r, L, R, v);
}
pushup(u);
}
void modify(int l, int r, const Tag &v) {
rangeModify(1, 0, n, l, r, v);
}
Info rangeQuery(int u, int l, int r, int L, int R) {
if (l >= L && r <= R) {
return info[u];
}
pushdown(u);
int mid = l + r >> 1;
Info x = Info(), y = Info();
if (mid >= L) {
x = rangeQuery(u << 1, l, mid, L, R);
}
if (mid < R) {
y = rangeQuery(u << 1 | 1, mid + 1, r, L, R);
}
return x + y;
}
Info query(int x) {
return rangeQuery(1, 0, n, x, x);
}
Info query(int l, int r) {
return rangeQuery(1, 0, n, l, r);
}
};
struct Tag {
int add;
void apply(Tag t) {
add += t.add;
}
};
struct Info {
int v;
void apply(Tag t) {
if (t.add == -1) {
v = 0;
}
}
};
Info operator+(const Info &a, const Info &b) {
Info c;
c.v = a.v + b.v;
return c;
}
Segment Tree 2
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
struct Node {
int l, r, val, tag;
};
struct Seg {
int n;
std::vector<Node> tr;
std::vector<int> a;
Seg(int _n) {
n = _n;
tr.resize(4 * n + 1);
a.resize(n + 1);
}
void pushup(int rt) {
pushdown(rt * 2), pushdown(rt * 2 + 1);
tr[rt].val = std::min(tr[rt * 2].val, tr[rt * 2 + 1].val);
}
void pushdown(int rt) {
if (tr[rt].tag == 0) return;
tr[rt].val += tr[rt].tag;
if (tr[rt].l != tr[rt].r) {
tr[rt * 2].tag += tr[rt].tag;
tr[rt * 2 + 1].tag += tr[rt].tag;
}
tr[rt].tag = 0;
}
void build(int rt, int l, int r) {
tr[rt].l = l, tr[rt].r = r;
if (l == r) {
tr[rt].val = a[l];
return;
}
int mid = (l + r) >> 1;
build(rt * 2, l, mid);
build(rt * 2 + 1, mid + 1, r);
pushup(rt);
}
void update(int rt, int l, int r, int x) {
pushdown(rt);
if (tr[rt].l == l && r == tr[rt].r) {
tr[rt].tag = x;
return;
}
int mid = (tr[rt].l + tr[rt].r) >> 1;
if (r <= mid) update(rt * 2, l, r, x);
else if (l > mid) update(rt * 2 + 1, l, r, x);
else update(rt * 2, l, mid, x), update(rt * 2 + 1, mid + 1, r, x);
pushup(rt);
}
int query(int rt, int l, int r) {
pushdown(rt);
if (tr[rt].l == l && tr[rt].r == r) {
return tr[rt].val;
}
int mid = (tr[rt].l + tr[rt].r) >> 1;
if (r <= mid) return query(rt * 2, l, r);
else if (l > mid) return query(rt * 2 + 1, l, r);
else return std::min(query(rt * 2, l, mid), query(rt * 2 + 1, mid + 1, r));
}
};
Graph Theory
DSU
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
struct DSU {
std::vector<int> f, siz;
DSU() {}
DSU(int n) {
init(n);
}
void init(int n) {
f.resize(n);
std::iota(f.begin(), f.end(), 0);
siz.assign(n, 1);
}
int find(int x) {
while (x != f[x]) {
x = f[x] = f[f[x]];
}
return x;
}
bool same(int x, int y) {
return find(x) == find(y);
}
bool merge(int x, int y) {
x = find(x);
y = find(y);
if (x == y) {
return false;
}
siz[x] += siz[y];
f[y] = x;
return true;
}
int size(int x) {
return siz[find(x)];
}
};
DSU with roll back
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
struct DSU {
std::vector<int> siz;
std::vector<int> f;
std::vector<std::array<int, 2>> his;
DSU(int n) : siz(n + 1, 1), f(n + 1) {
std::iota(f.begin(), f.end(), 0);
}
int find(int x) {
while (f[x] != x) {
x = f[x];
}
return x;
}
bool merge(int x, int y) {
x = find(x);
y = find(y);
if (x == y) {
return false;
}
if (siz[x] < siz[y]) {
std::swap(x, y);
}
his.push_back({x, y});
siz[x] += siz[y];
f[y] = x;
return true;
}
int time() {
return his.size();
}
void revert(int tm) {
while (his.size() > tm) {
auto [x, y] = his.back();
his.pop_back();
f[y] = y;
siz[x] -= siz[y];
}
}
};
Network Flow
ISAP
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
const int N=100;//点数，从1开始编号
const int M=100000;//边数
const int inf=1000000000;//无穷大
struct E{int t,f;E*nxt,*pair;}*g[N],*d[N],pool[M*2],*cur=pool;
int n,m,i,S,T,h[N],gap[N],maxflow;
void add(int s,int t,int f){//添加s->t的一条边，容量为f
E*p=cur++;p->t=t;p->f=f;p->nxt=g[s];g[s]=p;
p=cur++;p->t=s;p->f=0;p->nxt=g[t];g[t]=p;
g[s]->pair=g[t];g[t]->pair=g[s];
}
int sap(int v,int flow){
if(v==T)return flow;
int rec=0;
for(E*p=d[v];p;p=p->nxt)if(h[v]==h[p->t]+1&&p->f){
int ret=sap(p->t,min(flow-rec,p->f));
p->f-=ret;p->pair->f+=ret;d[v]=p;
if((rec+=ret)==flow)return flow;
}
if(!(--gap[h[v]]))h[S]=T;
gap[++h[v]]++;d[v]=g[v];
return rec;
}
int main(){
S=n+1;//源点
T=S+1;//汇点
for(cur=pool,i=1;i<=T;i++)g[i]=d[i]=NULL,h[i]=gap[i]=0;//初始化
//请在这里连边
for(gap[maxflow=0]=T,i=1;i<=T;i++)d[i]=g[i];
while(h[S]<T)maxflow+=sap(S,inf);
printf("maxflow=%d\n",maxflow);
}
MaxFlow
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
constexpr int inf = 1E9;
template<class T>
struct MaxFlow {
struct _Edge {
int to;
T cap;
_Edge(int to, T cap) : to(to), cap(cap) {}
};
int n;
std::vector<_Edge> e;
std::vector<std::vector<int>> g;
std::vector<int> cur, h;
MaxFlow() {}
MaxFlow(int n) {
init(n);
}
void init(int n) {
this->n = n;
e.clear();
g.assign(n, {});
cur.resize(n);
h.resize(n);
}
bool bfs(int s, int t) {
h.assign(n, -1);
std::queue<int> que;
h[s] = 0;
que.push(s);
while (!que.empty()) {
const int u = que.front();
que.pop();
for (int i : g[u]) {
auto [v, c] = e[i];
if (c > 0 && h[v] == -1) {
h[v] = h[u] + 1;
if (v == t) {
return true;
}
que.push(v);
}
}
}
return false;
}
T dfs(int u, int t, T f) {
if (u == t) {
return f;
}
auto r = f;
for (int &i = cur[u]; i < int(g[u].size()); ++i) {
const int j = g[u][i];
auto [v, c] = e[j];
if (c > 0 && h[v] == h[u] + 1) {
auto a = dfs(v, t, std::min(r, c));
e[j].cap -= a;
e[j ^ 1].cap += a;
r -= a;
if (r == 0) {
return f;
}
}
}
return f - r;
}
void addEdge(int u, int v, T c) {
g[u].push_back(e.size());
e.emplace_back(v, c);
g[v].push_back(e.size());
e.emplace_back(u, 0);
}
T flow(int s, int t) {
T ans = 0;
while (bfs(s, t)) {
cur.assign(n, 0);
ans += dfs(s, t, std::numeric_limits<T>::max());
}
return ans;
}
std::vector<bool> minCut() {
std::vector<bool> c(n);
for (int i = 0; i < n; i++) {
c[i] = (h[i] != -1);
}
return c;
}
struct Edge {
int from;
int to;
T cap;
T flow;
};
std::vector<Edge> edges() {
std::vector<Edge> a;
for (int i = 0; i < e.size(); i += 2) {
Edge x;
x.from = e[i + 1].to;
x.to = e[i].to;
x.cap = e[i].cap + e[i + 1].cap;
x.flow = e[i + 1].cap;
a.push_back(x);
}
return a;
}
};
MinCostFlow
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
template<class T>
struct MinCostFlow {
struct _Edge {
int to;
T cap;
T cost;
_Edge(int to_, T cap_, T cost_) : to(to_), cap(cap_), cost(cost_) {}
};
int n;
std::vector<_Edge> e;
std::vector<std::vector<int>> g;
std::vector<T> h, dis;
std::vector<int> pre;
bool dijsktra(int s, int t) {
dis.assign(n, std::numeric_limits<T>::max());
pre.assign(n, -1);
std::priority_queue<std::pair<T, int>, std::vector<std::pair<T, int>>, std::greater<std::pair<T, int>>> que;
dis[s] = 0;
que.emplace(0, s);
while (!que.empty()) {
T d = que.top().first;
int u = que.top().second;
que.pop();
if (dis[u] != d) continue;
for (int i : g[u]) {
int v = e[i].to;
T cap = e[i].cap;
T cost = e[i].cost;
if (cap > 0 && dis[v] > d + h[u] - h[v] + cost) {
dis[v] = d + h[u] - h[v] + cost;
pre[v] = i;
que.emplace(dis[v], v);
}
}
}
return dis[t] != std::numeric_limits<T>::max();
}
MinCostFlow() {}
MinCostFlow(int n_) {
init(n_);
}
void init(int n_) {
n = n_;
e.clear();
g.assign(n, {});
}
void addEdge(int u, int v, T cap, T cost) {
g[u].push_back(e.size());
e.emplace_back(v, cap, cost);
g[v].push_back(e.size());
e.emplace_back(u, 0, -cost);
}
std::pair<T, T> flow(int s, int t) {
T flow = 0;
T cost = 0;
h.assign(n, 0);
while (dijsktra(s, t)) {
for (int i = 0; i < n; ++i) {
h[i] += dis[i];
}
T aug = std::numeric_limits<int>::max();
for (int i = t; i != s; i = e[pre[i] ^ 1].to) {
aug = std::min(aug, e[pre[i]].cap);
}
for (int i = t; i != s; i = e[pre[i] ^ 1].to) {
e[pre[i]].cap -= aug;
e[pre[i] ^ 1].cap += aug;
}
flow += aug;
cost += aug * h[t];
}
return std::make_pair(flow, cost);
}
struct Edge {
int from;
int to;
T cap;
T cost;
T flow;
};
std::vector<Edge> edges() {
std::vector<Edge> a;
for (int i = 0; i < e.size(); i += 2) {
Edge x;
x.from = e[i + 1].to;
x.to = e[i].to;
x.cap = e[i].cap + e[i + 1].cap;
x.cost = e[i].cost;
x.flow = e[i + 1].cap;
a.push_back(x);
}
return a;
}
};
HLD(树链剖分)
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
struct HLD {
int n;
std::vector<int> siz, top, dep, parent, in, out, seq;
std::vector<std::vector<int>> adj;
int cur;
HLD() {}
HLD(int n) {
init(n);
}
void init(int n) {
this->n = n;
siz.resize(n);
top.resize(n);
dep.resize(n);
parent.resize(n);
in.resize(n);
out.resize(n);
seq.resize(n);
cur = 0;
adj.assign(n, {});
}
void addEdge(int u, int v) {
adj[u].push_back(v);
adj[v].push_back(u);
}
void work(int root = 0) {
top[root] = root;
dep[root] = 0;
parent[root] = -1;
dfs1(root);
dfs2(root);
}
void dfs1(int u) {
if (parent[u] != -1) {
adj[u].erase(std::find(adj[u].begin(), adj[u].end(), parent[u]));
}
siz[u] = 1;
for (auto &v : adj[u]) {
parent[v] = u;
dep[v] = dep[u] + 1;
dfs1(v);
siz[u] += siz[v];
if (siz[v] > siz[adj[u][0]]) {
std::swap(v, adj[u][0]);
}
}
}
void dfs2(int u) {
in[u] = cur++;
seq[in[u]] = u;
for (auto v : adj[u]) {
top[v] = v == adj[u][0] ? top[u] : v;
dfs2(v);
}
out[u] = cur;
}
int lca(int u, int v) {
while (top[u] != top[v]) {
if (dep[top[u]] > dep[top[v]]) {
u = parent[top[u]];
} else {
v = parent[top[v]];
}
}
return dep[u] < dep[v] ? u : v;
}
int dist(int u, int v) {
return dep[u] + dep[v] - 2 * dep[lca(u, v)];
}
int jump(int u, int k) {
if (dep[u] < k) {
return -1;
}
int d = dep[u] - k;
while (dep[top[u]] > d) {
u = parent[top[u]];
}
return seq[in[u] - dep[u] + d];
}
bool isAncester(int u, int v) {
return in[u] <= in[v] && in[v] < out[u];
}
int rootedParent(int u, int v) {
std::swap(u, v);
if (u == v) {
return u;
}
if (!isAncester(u, v)) {
return parent[u];
}
auto it = std::upper_bound(adj[u].begin(), adj[u].end(), v, [&](int x, int y) {
return in[x] < in[y];
}) - 1;
return *it;
}
int rootedSize(int u, int v) {
if (u == v) {
return n;
}
if (!isAncester(v, u)) {
return siz[v];
}
return n - siz[rootedParent(u, v)];
}
int rootedLca(int a, int b, int c) {
return lca(a, b) ^ lca(b, c) ^ lca(c, a);
}
};
SCC(强连通分量)
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
struct SCC {
static constexpr int N = 100005;
static constexpr int M = 200005;
struct Edge {
int u, v, nex;
}gra[M], graInv[M];
int n, m, ptr, ptrInv, ptrQue, head[N], headInv[N], que[N], f[N];
bitset<N> vis;
void init(int n, int m) {
this->n = n;
this->m = m;
}
void addEdge(int u, int v) {
++ptr;
gra[ptr].u = u;
gra[ptr].v = v;
gra[ptr].nex = head[u];
head[u] = ptr;
++ptrInv;
graInv[ptrInv].u = v;
graInv[ptrInv].v = u;
graInv[ptrInv].nex = headInv[v];
headInv[v] = ptrInv;
}
void dfs1(int x) {
vis[x] = 1;
for (int i = headInv[x]; i; i = graInv[i].nex) {
int v = graInv[i].v;
if (!vis.test(v)) dfs1(v);
}
que[++ptrQue] = x;
}
void dfs2(int x, int y) {
vis[x] = 0; f[x] = y;
res[y] = max(res[y], x);
for (int i = head[x]; i; i = gra[i].nex) {
int v = gra[i].v;
if (vis.test(v)) {
dfs2(v, y);
}
}
}
void getSCC() {
for (int i = 1; i <= n; ++i) {
if (!vis.test(i)) dfs1(i);
}
for (int i = n; i >= 1; --i) {
if (vis.test(que[i])) dfs2(que[i], que[i]);
}
}
void rebuild() {
ptr = 0;
memset(head, 0, sizeof(head));
for (int i = 1; i <= m; ++i) {
int u = gra[i].u, v = gra[i].v;
if (f[u] == f[v]) continue;
addEdge(f[u], f[v]);
}
}
};
Math
NTT
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
namespace NTT
{
const long long g = 3;
const long long p = 998244353;
long long wn[35];
long long pow2(long long a, long long b)
{
long long res = 1;
while (b)
{
if (b & 1) res = res * a % p;
a = a * a % p;
b >>= 1;
}
return res;
}
void getwn()
{
for (int i = 0; i < 25; i++) wn[i] = pow2(g, (p - 1) / (1LL << i));
}
void ntt(long long *a, int len, int f)
{
long long i, j = 0, t, k, w, id;
for (i = 1; i < len - 1; i++)
{
for (t = len; j ^= t >>= 1, ~j & t;);
if (i < j) swap(a[i], a[j]);
}
for (i = 1, id = 1; i < len; i <<= 1, id++)
{
t = i << 1;
for (j = 0; j < len; j += t)
{
for (k = 0, w = 1; k < i; k++, w = w * wn[id] % p)
{
long long x = a[j + k], y = w * a[j + k + i] % p;
a[j + k] = (x + y) % p;
a[j + k + i] = (x - y + p) % p;
}
}
}
if (f)
{
for (i = 1, j = len - 1; i < j; i++, j--) swap(a[i], a[j]);
long long inv = pow2(len, p - 2);
for (i = 0; i < len; i++) a[i] = a[i] * inv % p;
}
}
void mul(long long *a, long long *b, int l1, int l2)
{
int len, i;
for (len = 1; len <= l1 + l2; len <<= 1);
for (i = l1 + 1; i <= len; i++) a[i] = 0;
for (i = l2 + 1; i <= len; i++) b[i] = 0;
ntt(a, len, 0); ntt(b, len, 0);
for (i = 0; i < len; i++) a[i] = a[i] * b[i] % p;
ntt(a, len, 1);
}
};
GCD
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
template<typename T>
static T gcd(T a, T b) {
return b == 0 ? a : gcd(b, a % b);
}
template<typename T>
static T exgcd(T a, T b, T& x, T& y) {
if (b == 0) {
x = 1, y = 0;
return a;
}
T gcd = exgcd(b, a % b, y, x);
y -= a / b * x;
return gcd;
}
template<typename T>
static T MinSol(T a, T b, T c) {
T x, y, z;
z = exgcd(a, b, x, y);
if (c % z != 0) return -1;
x *= c / z;
b /= z;
if (b < 0) b = -b;
return (x % b + b) % b;
}
ModInt
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
template<class T>
constexpr T power(T a, u64 b, T res = 1) {
for (; b != 0; b /= 2, a *= a) {
if (b & 1) {
res *= a;
}
}
return res;
}
template<u32 P>
constexpr u32 mulMod(u32 a, u32 b) {
return u64(a) * b % P;
}
template<u64 P>
constexpr u64 mulMod(u64 a, u64 b) {
u64 res = a * b - u64(1.L * a * b / P - 0.5L) * P;
res %= P;
return res;
}
constexpr i64 safeMod(i64 x, i64 m) {
x %= m;
if (x < 0) {
x += m;
}
return x;
}
constexpr std::pair<i64, i64> invGcd(i64 a, i64 b) {
a = safeMod(a, b);
if (a == 0) {
return {b, 0};
}
i64 s = b, t = a;
i64 m0 = 0, m1 = 1;
while (t) {
i64 u = s / t;
s -= t * u;
m0 -= m1 * u;
std::swap(s, t);
std::swap(m0, m1);
}
if (m0 < 0) {
m0 += b / s;
}
return {s, m0};
}
template<std::unsigned_integral U, U P>
struct ModIntBase {
public:
constexpr ModIntBase() : x(0) {}
template<std::unsigned_integral T>
constexpr ModIntBase(T x_) : x(x_ % mod()) {}
template<std::signed_integral T>
constexpr ModIntBase(T x_) {
using S = std::make_signed_t<U>;
S v = x_ % S(mod());
if (v < 0) {
v += mod();
}
x = v;
}
constexpr static U mod() {
return P;
}
constexpr U val() const {
return x;
}
constexpr ModIntBase operator-() const {
ModIntBase res;
res.x = (x == 0 ? 0 : mod() - x);
return res;
}
constexpr ModIntBase inv() const {
return power(*this, mod() - 2);
}
constexpr ModIntBase &operator*=(const ModIntBase &rhs) & {
x = mulMod<mod()>(x, rhs.val());
return *this;
}
constexpr ModIntBase &operator+=(const ModIntBase &rhs) & {
x += rhs.val();
if (x >= mod()) {
x -= mod();
}
return *this;
}
constexpr ModIntBase &operator-=(const ModIntBase &rhs) & {
x -= rhs.val();
if (x >= mod()) {
x += mod();
}
return *this;
}
constexpr ModIntBase &operator/=(const ModIntBase &rhs) & {
return *this *= rhs.inv();
}
friend constexpr ModIntBase operator*(ModIntBase lhs, const ModIntBase &rhs) {
lhs *= rhs;
return lhs;
}
friend constexpr ModIntBase operator+(ModIntBase lhs, const ModIntBase &rhs) {
lhs += rhs;
return lhs;
}
friend constexpr ModIntBase operator-(ModIntBase lhs, const ModIntBase &rhs) {
lhs -= rhs;
return lhs;
}
friend constexpr ModIntBase operator/(ModIntBase lhs, const ModIntBase &rhs) {
lhs /= rhs;
return lhs;
}
friend constexpr std::istream &operator>>(std::istream &is, ModIntBase &a) {
i64 i;
is >> i;
a = i;
return is;
}
friend constexpr std::ostream &operator<<(std::ostream &os, const ModIntBase &a) {
return os << a.val();
}
friend constexpr bool operator==(const ModIntBase &lhs, const ModIntBase &rhs) {
return lhs.val() == rhs.val();
}
friend constexpr std::strong_ordering operator<=>(const ModIntBase &lhs, const ModIntBase &rhs) {
return lhs.val() <=> rhs.val();
}
private:
U x;
};
template<u32 P>
using ModInt = ModIntBase<u32, P>;
template<u64 P>
using ModInt64 = ModIntBase<u64, P>;
struct Barrett {
public:
Barrett(u32 m_) : m(m_), im((u64)(-1) / m_ + 1) {}
constexpr u32 mod() const {
return m;
}
constexpr u32 mul(u32 a, u32 b) const {
u64 z = a;
z *= b;
u64 x = u64((u128(z) * im) >> 64);
u32 v = u32(z - x * m);
if (m <= v) {
v += m;
}
return v;
}
private:
u32 m;
u64 im;
};
template<u32 Id>
struct DynModInt {
public:
constexpr DynModInt() : x(0) {}
template<std::unsigned_integral T>
constexpr DynModInt(T x_) : x(x_ % mod()) {}
template<std::signed_integral T>
constexpr DynModInt(T x_) {
int v = x_ % int(mod());
if (v < 0) {
v += mod();
}
x = v;
}
constexpr static void setMod(u32 m) {
bt = m;
}
static u32 mod() {
return bt.mod();
}
constexpr u32 val() const {
return x;
}
constexpr DynModInt operator-() const {
DynModInt res;
res.x = (x == 0 ? 0 : mod() - x);
return res;
}
constexpr DynModInt inv() const {
auto v = invGcd(x, mod());
assert(v.first == 1);
return v.second;
}
constexpr DynModInt &operator*=(const DynModInt &rhs) & {
x = bt.mul(x, rhs.val());
return *this;
}
constexpr DynModInt &operator+=(const DynModInt &rhs) & {
x += rhs.val();
if (x >= mod()) {
x -= mod();
}
return *this;
}
constexpr DynModInt &operator-=(const DynModInt &rhs) & {
x -= rhs.val();
if (x >= mod()) {
x += mod();
}
return *this;
}
constexpr DynModInt &operator/=(const DynModInt &rhs) & {
return *this *= rhs.inv();
}
friend constexpr DynModInt operator*(DynModInt lhs, const DynModInt &rhs) {
lhs *= rhs;
return lhs;
}
friend constexpr DynModInt operator+(DynModInt lhs, const DynModInt &rhs) {
lhs += rhs;
return lhs;
}
friend constexpr DynModInt operator-(DynModInt lhs, const DynModInt &rhs) {
lhs -= rhs;
return lhs;
}
friend constexpr DynModInt operator/(DynModInt lhs, const DynModInt &rhs) {
lhs /= rhs;
return lhs;
}
friend constexpr std::istream &operator>>(std::istream &is, DynModInt &a) {
i64 i;
is >> i;
a = i;
return is;
}
friend constexpr std::ostream &operator<<(std::ostream &os, const DynModInt &a) {
return os << a.val();
}
friend constexpr bool operator==(const DynModInt &lhs, const DynModInt &rhs) {
return lhs.val() == rhs.val();
}
friend constexpr std::strong_ordering operator<=>(const DynModInt &lhs, const DynModInt &rhs) {
return lhs.val() <=> rhs.val();
}
private:
u32 x;
static Barrett bt;
};
template<u32 Id>
Barrett DynModInt<Id>::bt = 998244353;
using Z = ModInt<998244353>;
MInt
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
template<class T>
constexpr T power(T a, i64 b) {
T res = 1;
while (b) {
if (b & 1) res *= a;
a *= a;
b >>= 1;
}
return res;
}
template<int P>
struct MInt {
int x;
constexpr MInt() : x{} {}
constexpr MInt(i64 x) : x{norm(x % getMod())} {}
static int Mod;
constexpr static int getMod() {
if (P > 0) {
return P;
} else {
return Mod;
}
}
constexpr static void setMod(int Mod_) {
Mod = Mod_;
}
constexpr int norm(int x) const {
if (x < 0) {
x += getMod();
}
if (x >= getMod()) {
x -= getMod();
}
return x;
}
constexpr int val() const {
return x;
}
explicit constexpr operator int() const {
return x;
}
constexpr MInt operator-() const {
MInt res;
res.x = norm(getMod() - x);
return res;
}
constexpr MInt inv() const {
assert(x != 0);
return power(*this, getMod() - 2);
}
constexpr MInt &operator*=(MInt rhs) & {
x = 1LL * x * rhs.x % getMod();
return *this;
}
constexpr MInt &operator+=(MInt rhs) & {
x = norm(x + rhs.x);
return *this;
}
constexpr MInt &operator-=(MInt rhs) & {
x = norm(x - rhs.x);
return *this;
}
constexpr MInt &operator/=(MInt rhs) & {
return *this *= rhs.inv();
}
friend constexpr MInt operator*(MInt lhs, MInt rhs) {
MInt res = lhs;
res *= rhs;
return res;
}
friend constexpr MInt operator+(MInt lhs, MInt rhs) {
MInt res = lhs;
res += rhs;
return res;
}
friend constexpr MInt operator-(MInt lhs, MInt rhs) {
MInt res = lhs;
res -= rhs;
return res;
}
friend constexpr MInt operator/(MInt lhs, MInt rhs) {
MInt res = lhs;
res /= rhs;
return res;
}
friend constexpr std::istream &operator>>(std::istream &is, MInt &a) {
i64 v;
is >> v;
a = MInt(v);
return is;
}
friend constexpr std::ostream &operator<<(std::ostream &os, const MInt &a) {
return os << a.val();
}
friend constexpr bool operator==(MInt lhs, MInt rhs) {
return lhs.val() == rhs.val();
}
friend constexpr bool operator!=(MInt lhs, MInt rhs) {
return lhs.val() != rhs.val();
}
friend constexpr bool operator<(MInt lhs, MInt rhs) {
return lhs.val() < rhs.val();
}
};
template<>
int MInt<0>::Mod = 1;
template<int V, int P>
constexpr MInt<P> CInv = MInt<P>(V).inv();
constexpr int P = 998244353;
using Z = MInt<P>;
MLong
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
template<class T>
constexpr T power(T a, i64 b) {
T res = 1;
for (; b; b /= 2, a *= a) {
if (b % 2) {
res *= a;
}
}
return res;
}
constexpr i64 mul(i64 a, i64 b, i64 p) {
i64 res = a * b - i64(1.L * a * b / p) * p;
res %= p;
if (res < 0) {
res += p;
}
return res;
}
template<i64 P>
struct MLong {
i64 x;
constexpr MLong() : x{} {}
constexpr MLong(i64 x) : x{norm(x % getMod())} {}
static i64 Mod;
constexpr static i64 getMod() {
if (P > 0) {
return P;
} else {
return Mod;
}
}
constexpr static void setMod(i64 Mod_) {
Mod = Mod_;
}
constexpr i64 norm(i64 x) const {
if (x < 0) {
x += getMod();
}
if (x >= getMod()) {
x -= getMod();
}
return x;
}
constexpr i64 val() const {
return x;
}
explicit constexpr operator i64() const {
return x;
}
constexpr MLong operator-() const {
MLong res;
res.x = norm(getMod() - x);
return res;
}
constexpr MLong inv() const {
assert(x != 0);
return power(*this, getMod() - 2);
}
constexpr MLong &operator*=(MLong rhs) & {
x = mul(x, rhs.x, getMod());
return *this;
}
constexpr MLong &operator+=(MLong rhs) & {
x = norm(x + rhs.x);
return *this;
}
constexpr MLong &operator-=(MLong rhs) & {
x = norm(x - rhs.x);
return *this;
}
constexpr MLong &operator/=(MLong rhs) & {
return *this *= rhs.inv();
}
friend constexpr MLong operator*(MLong lhs, MLong rhs) {
MLong res = lhs;
res *= rhs;
return res;
}
friend constexpr MLong operator+(MLong lhs, MLong rhs) {
MLong res = lhs;
res += rhs;
return res;
}
friend constexpr MLong operator-(MLong lhs, MLong rhs) {
MLong res = lhs;
res -= rhs;
return res;
}
friend constexpr MLong operator/(MLong lhs, MLong rhs) {
MLong res = lhs;
res /= rhs;
return res;
}
friend constexpr std::istream &operator>>(std::istream &is, MLong &a) {
i64 v;
is >> v;
a = MLong(v);
return is;
}
friend constexpr std::ostream &operator<<(std::ostream &os, const MLong &a) {
return os << a.val();
}
friend constexpr bool operator==(MLong lhs, MLong rhs) {
return lhs.val() == rhs.val();
}
friend constexpr bool operator!=(MLong lhs, MLong rhs) {
return lhs.val() != rhs.val();
}
};
template<>
i64 MLong<0LL>::Mod = i64(1E18) + 9;
template<int P>
struct MInt {
int x;
constexpr MInt() : x{} {}
constexpr MInt(i64 x) : x{norm(x % getMod())} {}
static int Mod;
constexpr static int getMod() {
if (P > 0) {
return P;
} else {
return Mod;
}
}
constexpr static void setMod(int Mod_) {
Mod = Mod_;
}
constexpr int norm(int x) const {
if (x < 0) {
x += getMod();
}
if (x >= getMod()) {
x -= getMod();
}
return x;
}
constexpr int val() const {
return x;
}
explicit constexpr operator int() const {
return x;
}
constexpr MInt operator-() const {
MInt res;
res.x = norm(getMod() - x);
return res;
}
constexpr MInt inv() const {
assert(x != 0);
return power(*this, getMod() - 2);
}
constexpr MInt &operator*=(MInt rhs) & {
x = 1LL * x * rhs.x % getMod();
return *this;
}
constexpr MInt &operator+=(MInt rhs) & {
x = norm(x + rhs.x);
return *this;
}
constexpr MInt &operator-=(MInt rhs) & {
x = norm(x - rhs.x);
return *this;
}
constexpr MInt &operator/=(MInt rhs) & {
return *this *= rhs.inv();
}
friend constexpr MInt operator*(MInt lhs, MInt rhs) {
MInt res = lhs;
res *= rhs;
return res;
}
friend constexpr MInt operator+(MInt lhs, MInt rhs) {
MInt res = lhs;
res += rhs;
return res;
}
friend constexpr MInt operator-(MInt lhs, MInt rhs) {
MInt res = lhs;
res -= rhs;
return res;
}
friend constexpr MInt operator/(MInt lhs, MInt rhs) {
MInt res = lhs;
res /= rhs;
return res;
}
friend constexpr std::istream &operator>>(std::istream &is, MInt &a) {
i64 v;
is >> v;
a = MInt(v);
return is;
}
friend constexpr std::ostream &operator<<(std::ostream &os, const MInt &a) {
return os << a.val();
}
friend constexpr bool operator==(MInt lhs, MInt rhs) {
return lhs.val() == rhs.val();
}
friend constexpr bool operator!=(MInt lhs, MInt rhs) {
return lhs.val() != rhs.val();
}
};
template<>
int MInt<0>::Mod = 998244353;
template<int V, int P>
constexpr MInt<P> CInv = MInt<P>(V).inv();
constexpr i64 P = i64(1E18) + 9;
constexpr i64 B = 114514;
using Z = MLong<P>;
Prime
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
std::vector<int> minp, primes;
void sieve(int n) {
minp.assign(n + 1, 0);
for (int i = 2; i <= n; ++i) {
if (minp[i] == 0) {
minp[i] = i;
primes.push_back(i);
}
for (auto p : primes) {
if (i * p > n)
break;
minp[i * p] = p;
if (p == minp[i])
break;
}
}
}
bool isPrime(int x) {
return x == minp[x];
}
QuickCal
1
2
3
4
5
6
7
8
9
10
11
12
template<class T>
struct QuickCal {
constexpr T power(T a, i64 b) {
T res = 1;
while (b) {
if (b & 1) res *= a;
a *= a;
b >>= 1;
}
return res;
}
};
Polynomial
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
vector<int> rev;
vector<Z> roots{0, 1};
void dft(std::vector<Z> &a) {
int n = a.size();
if (int(rev.size()) != n) {
int k = __builtin_ctz(n) - 1;
rev.resize(n);
for (int i = 0; i < n; i++) {
rev[i] = rev[i >> 1] >> 1 | (i & 1) << k;
}
}
for (int i = 0; i < n; i++) {
if (rev[i] < i) {
std::swap(a[i], a[rev[i]]);
}
}
if (int(roots.size()) < n) {
int k = __builtin_ctz(roots.size());
roots.resize(n);
while ((1 << k) < n) {
Z e = power(Z(3), (P - 1) >> (k + 1));
for (int i = 1 << (k - 1); i < (1 << k); i++) {
roots[2 * i] = roots[i];
roots[2 * i + 1] = roots[i] * e;
}
k++;
}
}
for (int k = 1; k < n; k *= 2) {
for (int i = 0; i < n; i += 2 * k) {
for (int j = 0; j < k; j++) {
Z u = a[i + j];
Z v = a[i + j + k] * roots[k + j];
a[i + j] = u + v;
a[i + j + k] = u - v;
}
}
}
}
void idft(std::vector<Z> &a) {
int n = a.size();
std::reverse(a.begin() + 1, a.end());
dft(a);
Z inv = (1 - P) / n;
for (int i = 0; i < n; i++) {
a[i] *= inv;
}
}
struct Poly {
std::vector<Z> a;
Poly() {}
Poly(int x) : a(x) {}
Poly(const vector<Z> &_a) : a(_a) {}
Poly(const initializer_list<Z> &_a) : a(_a) {}
int size() const { return a.size();}
void resize(int n) { a.resize(n);}
Z operator[](int idx) const {
if (idx < size()) {
return a[idx];
} else {
return 0;
}
}
Z &operator[](int idx) {
return a[idx];
}
Poly mulxk(int k) const {
auto b = a;
b.insert(b.begin(), k, 0);
return Poly(b);
}
Poly modxk(int k) const {
k = std::min(k, size());
return Poly(std::vector<Z>(a.begin(), a.begin() + k));
}
Poly divxk(int k) const {
if (size() <= k) {
return Poly();
}
return Poly(std::vector<Z>(a.begin() + k, a.end()));
}
friend Poly operator+(const Poly &a, const Poly &b) {
std::vector<Z> res(std::max(a.size(), b.size()));
for (int i = 0; i < int(res.size()); i++) {
res[i] = a[i] + b[i];
}
return Poly(res);
}
friend Poly operator-(const Poly &a, const Poly &b) {
std::vector<Z> res(std::max(a.size(), b.size()));
for (int i = 0; i < int(res.size()); i++) {
res[i] = a[i] - b[i];
}
return Poly(res);
}
friend Poly operator*(Poly a, Poly b) {
if (a.size() == 0 || b.size() == 0) {
return Poly();
}
int sz = 1, tot = a.size() + b.size() - 1;
while (sz < tot) {
sz *= 2;
}
a.a.resize(sz);
b.a.resize(sz);
dft(a.a);
dft(b.a);
for (int i = 0; i < sz; ++i) {
a.a[i] = a[i] * b[i];
}
idft(a.a);
a.resize(tot);
return a;
}
friend Poly operator*(Z a, Poly b) {
for (int i = 0; i < int(b.size()); i++) {
b[i] *= a;
}
return b;
}
friend Poly operator*(Poly a, Z b) {
for (int i = 0; i < int(a.size()); i++) {
a[i] *= b;
}
return a;
}
Poly &operator+=(Poly b) {
return (*this) = (*this) + b;
}
Poly &operator-=(Poly b) {
return (*this) = (*this) - b;
}
Poly &operator*=(Poly b) {
return (*this) = (*this) * b;
}
Poly deriv() const {
if (a.empty()) {
return Poly();
}
std::vector<Z> res(size() - 1);
for (int i = 0; i < size() - 1; ++i) {
res[i] = (i + 1) * a[i + 1];
}
return Poly(res);
}
Poly integr() const {
std::vector<Z> res(size() + 1);
for (int i = 0; i < size(); ++i) {
res[i + 1] = a[i] / (i + 1);
}
return Poly(res);
}
Poly inv(int m) const {
Poly x{a[0].inv()};
int k = 1;
while (k < m) {
k *= 2;
x = (x * (Poly{2} - modxk(k) * x)).modxk(k);
}
return x.modxk(m);
}
//需要保证首项为1
Poly log(int m) const {
return (deriv() * inv(m)).integr().modxk(m);
}
//需要保证首项为0
Poly exp(int m) const {
Poly x{1};
int k = 1;
while (k < m) {
k *= 2;
x = (x * (Poly{1} - x.log(k) + modxk(k))).modxk(k);
}
return x.modxk(m);
}
Poly pow(int k, int m) const {
int i = 0;
while (i < size() && a[i].val() == 0) {
i++;
}
if (i == size() || 1LL * i * k >= m) {
return Poly(std::vector<Z>(m));
}
Z v = a[i];
auto f = divxk(i) * v.inv();
return (f.log(m - i * k) * Z(k)).exp(m - i * k).mulxk(i * k) * power(v, k);
}
Poly sqrt(int m) const {
Poly x{1};
int k = 1;
while (k < m) {
k *= 2;
x = (x + (modxk(k) * x.inv(k)).modxk(k)) * Z(((P + 1) / 2));
}
return x.modxk(m);
}
//转置卷积
Poly mulT(Poly b) const {
if (b.size() == 0) {
return Poly();
}
int n = b.size();
std::reverse(b.a.begin(), b.a.end());
return ((*this) * b).divxk(n - 1);
}
//多点求值，https://jkloverdcoi.github.io/2020/08/04/转置原理及其应用/
std::vector<Z> eval(std::vector<Z> x) const {
if (size() == 0) {
return std::vector<Z>(x.size(), 0);
}
const int n = std::max(int(x.size()), size());
std::vector<Poly> q(4 * n);
std::vector<Z> ans(x.size());
x.resize(n);
std::function<void(int, int, int)> build = [&](int p, int l, int r) {
if (r - l == 1) {
q[p] = Poly{1, -x[l]};
} else {
int m = (l + r) / 2;
build(2 * p, l, m);
build(2 * p + 1, m, r);
q[p] = q[2 * p] * q[2 * p + 1];
}
};
build(1, 0, n);
std::function<void(int, int, int, const Poly &)> work = [&](int p, int l, int r, const Poly &num) {
if (r - l == 1) {
if (l < int(ans.size())) {
ans[l] = num[0];
}
} else {
int m = (l + r) / 2;
work(2 * p, l, m, num.mulT(q[2 * p + 1]).modxk(m - l));
work(2 * p + 1, m, r, num.mulT(q[2 * p]).modxk(r - m));
}
};
work(1, 0, n, mulT(q[1].inv(n)));
return ans;
}
};
comb
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
struct Comb {
int n;
std::vector<Z> _fac;
std::vector<Z> _invfac;
std::vector<Z> _inv;
Comb() : n{0}, _fac{1}, _invfac{1}, _inv{0} {}
Comb(int n) : Comb() {
init(n);
}
void init(int m) {
m = std::min(m, Z::getMod() - 1);
if (m <= n) return;
_fac.resize(m + 1);
_invfac.resize(m + 1);
_inv.resize(m + 1);
for (int i = n + 1; i <= m; i++) {
_fac[i] = _fac[i - 1] * i;
}
_invfac[m] = _fac[m].inv();
for (int i = m; i > n; i--) {
_invfac[i - 1] = _invfac[i] * i;
_inv[i] = _invfac[i] * _fac[i - 1];
}
n = m;
}
Z fac(int m) {
if (m > n) init(2 * m);
return _fac[m];
}
Z invfac(int m) {
if (m > n) init(2 * m);
return _invfac[m];
}
Z inv(int m) {
if (m > n) init(2 * m);
return _inv[m];
}
Z binom(int n, int m) {
if (n < m || m < 0) return 0;
return fac(n) * invfac(m) * invfac(n - m);
}
} comb;
Euler
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
int phi(int n) {
int res = n;
for (int i = 2; i * i <= n; i++) {
if (n % i == 0) {
while (n % i == 0) {
n /= i;
}
res = res / i * (i - 1);
}
}
if (n > 1) {
res = res / n * (n - 1);
}
return res;
}
Other
random_int
1
2
3
4
5
6
int random_int(int l, int h) {
std::random_device rd;
std::mt19937 gen(rd());
std::uniform_int_distribution<> distrib(l, h);
return distrib(gen);
}
__int128
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
std::istream& operator>>(std::istream& in, __int128& x){
x = 0;
int w = 0;
char ch = 0;
while (!isdigit(ch)) {
w |= ch == '-';
ch = getchar();
}
while (isdigit(ch)) {
x = (x << 3) + (x << 1) + (ch ^ 48);
ch = getchar();
}
if (w) {
x = -x;
}
return in;
}
std::ostream& operator<<(std::ostream& os, __int128 x) {
if (!x) os << 0;
else if (x < 0) os << '-', x = -x;
std::string s;
while (x) s.push_back(x % 10 + '0'), x /= 10;
reverse(s.begin(), s.end());
os << s;
return os;
}
Fast Read
1
2
3
4
5
6
7
8
9
10
inline char nc() {
static char buf[100000], *p1 = buf,*p2 = buf;
return p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 100000, stdin), p1 == p2) ? EOF : *p1++;
}
inline int read(){
char ch = nc(); int sum = 0;
while(!(ch >= '0' && ch <= '9'))ch = nc();
while(ch >= '0' && ch <= '9') sum = sum * 10 + ch - 48,ch = nc();
return sum;
}
String
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
auto getPre = [](const std::string& s) {
int n = s.length();
std::vector<int> pre(n);
for (int i = 1; i < n; ++i) {
int j = pre[i - 1];
while (j && s[j] != s[i]) j = pre[j - 1];
if (s[i] == s[j]) pre[i] = j + 1;
}
return pre;
};
auto kmp = [&](const std::string& s1, const std::string& s2) {
int n1 = s1.length(), n2 = s2.length();
assert(n1 > 0 && n2 > 0);
std::vector<int> pre = getPre(s2);
std::vector<int> ans;
int j = 0;
for (int i = 0; i < n1; i++) {
while (j && s2[j] != s1[i]) j = pre[j - 1];
j += (s2[j] == s1[i]);
if (j == n2) ans.push_back(i - n2 + 1);
}
return ans;
};
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
struct AhoCorasickAutomaton {
static constexpr int N = 1e6 + 10;
int tr[N][26];
int val[N], fail[N], idx;
int ed[N];
void insert(string s, int id) {
int u = 0;
for (int i = 0; i < s.size(); i++) {
if (!tr[u][s[i] - 'a']) {
tr[u][s[i] - 'a'] = ++ idx;
}
u = tr[u][s[i] - 'a'];
}
ed[u] = id;
}
void build() {
queue<int> q;
for (int i = 0; i < 26; i++) {
if (tr[0][i]) {
q.push(tr[0][i]);
}
}
while (q.size()) {
int u = q.front();
q.pop();
for (int i = 0; i < 26; i++) {
if (tr[u][i]) {
fail[tr[u][i]] = tr[fail[u]][i];
q.push(tr[u][i]);
} else {
tr[u][i] = tr[fail[u]][i];
}
}
}
}
void FailTree() {
vector<int> adj[idx + 1];
for (int i = 1; i <= idx; i++) {
adj[fail[i]].push_back(i);
}
function<void(int)> dfs = [&](int u) {
for (auto v : adj[u]) {
dfs(v);
val[u] += val[v];
}
};
dfs(0);
}
}AC;
SA
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
constexpr int N = 1e6 + 10;
char s[N];
int n, sa[N], rk[N], oldrk[N << 1];
int id[N], key1[N], cnt[N], height[N];
bool cmp(int x, int y, int w) {
return oldrk[x] == oldrk[y] && oldrk[x + w] == oldrk[y + w];
}
void getSA() {
int m = 127;
for (int i = 1; i <= n; i++) {
cnt[rk[i] = s[i]] ++;
}
for (int i = 1; i <= m; i++) {
cnt[i] += cnt[i - 1];
}
for (int i = n; i >= 1; i--) {
sa[cnt[rk[i]]--] = i;
}
for (int w = 1, i, p; w < n; w <<= 1, m = p) {
for (p = 0, i = n; i > n - w; i--) {
id[++ p] = i;
}
for (i = 1; i <= n; i++) {
if (sa[i] > w) {
id[++ p] = sa[i] - w;
}
}
memset(cnt, 0, sizeof cnt);
for (i = 1; i <= n; i++) {
++ cnt[key1[i] = rk[id[i]]];
}
for (i = 1; i <= m; i++) {
cnt[i] += cnt[i - 1];
}
for (i = n; i >= 1; i--) {
sa[cnt[key1[i]] --] = id[i];
}
memcpy(oldrk + 1, rk + 1, n * sizeof (int));
for (p = 0, i = 1; i <= n; i++) {
rk[sa[i]] = cmp(sa[i], sa[i - 1], w) ? p : ++ p;
}
}
}
void getHeight() {
for (int i = 1, k = 0; i <= n; i++) {
if (k) {
k --;
}
while (s[i + k] == s[sa[rk[i] - 1] + k]) {
k ++;
}
height[rk[i]] = k;
}
}
SAM
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
struct SuffixAutomaton {
static constexpr int ALPHABET_SIZE = 26, N = 1e5;
struct Node {
int len;
int link;
int next[ALPHABET_SIZE];
Node() : len(0), link(0), next{} {}
} t[2 * N];
int cntNodes;
SuffixAutomaton() {
cntNodes = 1;
fill(t[0].next, t[0].next + ALPHABET_SIZE, 1);
t[0].len = -1;
}
int extend(int p, int c) {
if (t[p].next[c]) {
int q = t[p].next[c];
if (t[q].len == t[p].len + 1) return q;
int r = ++cntNodes;
t[r].len = t[p].len + 1;
t[r].link = t[q].link;
copy(t[q].next, t[q].next + ALPHABET_SIZE, t[r].next);
t[q].link = r;
while (t[p].next[c] == q) {
t[p].next[c] = r;
p = t[p].link;
}
return r;
}
int cur = ++cntNodes;
t[cur].len = t[p].len + 1;
while (!t[p].next[c]) {
t[p].next[c] = cur;
p = t[p].link;
}
t[cur].link = extend(p, c);
return cur;
}
}SAM;
PAM
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
struct PalindromicAutomaton {
static constexpr int N = 5e5 + 10;
int cnt[N], fail[N], len[N], ch[N][26];
int idx, last, sz;
char s[N];
int node(int l) {
sz ++, len[sz] = l;
fail[sz] = cnt[sz] = 0;
return sz;
}
void init() {
sz = -1, last = 0;
node(0), node(-1);
memset(ch, 0, sizeof ch);
idx = 0;
s[idx] = '$';
fail[0] = 1;
}
int getfail(int x) {
while (s[idx - len[x] - 1] != s[idx]) {
x = fail[x];
}
return x;
}
void insert(char c) {
s[++ idx] = c;
int now = getfail(last);
if (!ch[now][c - 'a']) {
int x = node(len[now] + 2);
fail[x] = ch[getfail(fail[now])][c - 'a'];
ch[now][c - 'a'] = x;
}
last = ch[now][c - 'a'];
cnt[last] ++;
}
void solve() {
ll ans = 0;
for (int i = sz; i >= 0; i--) {
cnt[fail[i]] += cnt[i];
}
}
}pam;
Manacher
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
std::vector<int> manacher(std::string s) {
std::string t = "#";
for (auto c : s) {
t += c;
t += '#';
}
int n = t.size();
std::vector<int> r(n);
for (int i = 0, j = 0; i < n; i++) {
if (2 * j - i >= 0 && j + r[j] > i) {
r[i] = std::min(r[2 * j - i], j + r[j] - i);
}
while (i - r[i] >= 0 && i + r[i] < n && t[i - r[i]] == t[i + r[i]]) {
r[i] += 1;
}
if (i + r[i] > j + r[j]) {
j = i;
}
}
return r;
}
Z function
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
std::vector<int> Z(std::string s) {
int n = s.size();
std::vector<int> z(n + 1);
z[0] = n;
for (int i = 1, j = 1; i < n; i++) {
z[i] = std::max(0, std::min(j + z[j] - i, z[i - j]));
while (i + z[i] < n && s[z[i]] == s[i + z[i]]) {
z[i]++;
}
if (i + z[i] > j + z[j]) {
j = i;
}
}
return z;
}
本文由作者按照 CC BY 4.0 进行授权