[BZOJ 3512]DZY Loves Math IV
给定正整数\(n, m\)，求\(\sum_{i = 1}^n\sum_{j = 1}^m \varphi(ij)\)，膜\(10^9+7\)输出。
\(n\leq 10^5, m\leq 10^9\)。
[BZOJ 4916]神犇与蒟蒻
给定\(n\)，分别求出\(\mu(x^2)\)和\(\varphi(x^2)\)的前\(n\)项和。
\(1\leq n\leq 10^9\)。
[BZOJ 4025]二分图
只会做水题力……qwq（搞错大小写，自裁请）
每条边存在的时间事一个区间，因此考虑线段树分治……首先当然要写一个可撤销并查集。
然后每个点维护他到父亲的边的边权膜2，合并啥的和食物链就很类似力……然后判定已联通两点间边数的奇偶性就直接把两个点到根的值加起来就行了，因为LCA到根的部分会算两遍，对答案无影响。
代码：
/**************************************************************
Problem: 4025
User: danihao123
Language: C++
Result: Accepted
Time:9928 ms
Memory:34840 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
#include <vector>
#include <stack>
const int maxn = 100005;
const int maxno = maxn << 2;
const int maxm = 200005;
int n, m, T;
int p[maxn], siz[maxn], d[maxn];
void init_dsu() {
for(int i = 1; i <= n; i ++) {
p[i] = i; siz[i] = 1; d[i] = 0;
}
}
int get_fa(int x) {
if(p[x] == x) {
return x;
} else {
return get_fa(p[x]);
}
}
int get_d(int x) {
if(p[x] == x) {
return 0;
} else {
return (d[x] + get_d(p[x])) % 2;
}
}
typedef std::pair<int, int> upd;
upd link_set(int x, int y, int v) {
if(siz[x] > siz[y]) std::swap(x, y);
p[x] = y; siz[y] += siz[x]; d[x] = v;
return std::make_pair(x, y);
}
upd merge_set(int x, int y) {
int v = (get_d(x) + get_d(y) + 1) % 2;
return link_set(get_fa(x), get_fa(y), v);
}
void unlink_set(const upd &pr) {
int x = pr.first, y = pr.second;
p[x] = x; siz[y] -= siz[x]; d[x] = 0;
}
bool is_same(int x, int y) {
return (get_fa(x) == get_fa(y));
}
std::vector<upd> event[maxno];
std::stack<std::pair<int, upd> > S;
void update(int o, int L, int R, int ql, int qr, upd v) {
if(ql <= L && R <= qr) {
event[o].push_back(v);
} else {
int M = (L + R) / 2;
if(ql <= M) update(o << 1, L, M, ql, qr, v);
if(qr > M) update(o << 1 | 1, M + 1, R, ql, qr, v);
}
}
int val[maxno];
void solve(int o, int L, int R) {
for(int i = 0; i < event[o].size(); i ++) {
const upd &pr = event[o][i];
int x = pr.first, y = pr.second;
if(is_same(x, y)) {
if((get_d(x) + get_d(y)) % 2 == 0) {
val[o] = 0;
break;
}
} else {
S.push(std::make_pair(o, merge_set(x, y)));
}
}
if(L == R) {
if(val[o] != 0) val[o] = 1;
} else {
if(val[o] != 0) {
int M = (L + R) / 2;
solve(o << 1, L, M);
solve(o << 1 | 1, M + 1, R);
}
}
while((!S.empty()) && S.top().first >= o) {
upd v = S.top().second; S.pop();
unlink_set(v);
}
}
void print(int o, int L, int R) {
if(L == R) {
if(val[o]) {
puts("Yes");
} else {
puts("No");
}
} else {
if(val[o] != -1) {
val[o << 1] = val[o];
val[o << 1 | 1] = val[o];
}
int M = (L + R) / 2;
print(o << 1, L, M);
print(o << 1 | 1, M + 1, R);
}
}
int main() {
memset(val, -1, sizeof(val));
scanf("%d%d%d", &n, &m, &T);
for(int i = 1; i <= m; i ++) {
int u, v, s, t; scanf("%d%d%d%d", &u, &v, &s, &t);
s ++; if(s <= t) update(1, 1, T, s, t, std::make_pair(u, v));
}
init_dsu(); solve(1, 1, T); print(1, 1, T);
return 0;
}
[BZOJ 2654]tree
APIO的时候听了一下wqs二分，做了这题，然后现在才写题解……
wqs二分的思路大抵就是你要求必须选\(k\)个，那就二分每个操作的一个“额外代价”，然后进行没有限制的选择。然后最后选出来的个数事和你二分的代价正相关/反相关的。
这道题的话，就二分选择黑边的代价，然后跑一般的最小生成树（有相同边权时要选择黑边！）。当然我们会遇到一个问题，就是二分到\(x\)的时候选的比\(k\)多，到\(x + 1\)的时候又比\(k\)少了。这道题的处理方法，事考虑代价为\(x\)时，其实存在选\(k\)个的最优解（如果说大于\(x\)那就没了），因此我们钦点代价为\(x\)跑一遍限制只能选\(k\)条黑边的最短路即可。
代码：
/**************************************************************
Problem: 2654
User: danihao123
Language: C++
Result: Accepted
Time:5756 ms
Memory:5856 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
#include <vector>
const int maxn = 50005;
const int maxm = 100005;
struct Edge {
int u, v, d;
int typ;
bool operator <(const Edge &res) const {
if(d == res.d) {
return typ < res.typ;
} else {
return d < res.d;
}
}
};
Edge E[maxm];
int n, m, need;
int par[maxn], siz[maxn];
void init_dsu() {
for(int i = 1; i <= n; i ++) {
par[i] = i;
siz[i] = 1;
}
}
int get_fa(int x) {
if(par[x] == x) {
return x;
} else {
return (par[x] = get_fa(par[x]));
}
}
void link_set(int x, int y) {
if(siz[x] > siz[y]) std::swap(x, y);
siz[y] += siz[x];
par[x] = y;
}
void merge_set(int x, int y) {
return link_set(get_fa(x), get_fa(y));
}
bool is_same(int x, int y) {
return (get_fa(x) == get_fa(y));
}
typedef long long ll;
int ans = 0x7fffffff;
bool kruskal(int delta) {
std::vector<Edge> vec;
for(int i = 1; i <= m; i ++) {
Edge e = E[i];
if(e.typ == 0) {
e.d += delta;
}
vec.push_back(e);
}
std::sort(vec.begin(), vec.end());
int used = 0; ll tot = -(ll(delta)) * (ll(need));
init_dsu();
for(int i = 0; i < m; i ++) {
const Edge &e = vec[i];
int u = e.u, v = e.v;
if(!is_same(u, v)) {
if(used == need && e.typ == 0) continue;
merge_set(u, v); tot += e.d;
if(e.typ == 0) used ++;
}
}
if(used == need) {
ans = std::min(ans, (int(tot)));
return true;
} else {
return false;
}
}
int main() {
scanf("%d%d%d", &n, &m, &need);
for(int i = 1; i <= m; i ++) {
scanf("%d%d%d%d", &E[i].u, &E[i].v, &E[i].d, &E[i].typ);
E[i].u ++; E[i].v ++;
}
int L = -10000001, R = 10000001;
while(true) {
#ifdef LOCAL
printf("Range (%d, %d)\n", L, R);
fflush(stdout);
#endif
if(R - L <= 3) {
for(int i = R; i >= L; i --) {
if(kruskal(i)) {
break;
}
}
break;
}
int M = L + (R - L) / 2;
if(kruskal(M)) {
L = M;
} else {
R = M;
}
}
printf("%d\n", ans);
return 0;
}
[BZOJ 1004][HNOI2008]Cards
肉肉肉肉，，，
碰置换群啥的不是第一次了……这个题就是给你一堆置换（不要忘记还有一个幺元啊），然后限制颜色数量，求本质不同解个数。那么考虑使用Burnside引理，接下来考虑怎么计算每个置换的不动点数量，这个要求每个循环的颜色一致（不就事Polya定理了吗），所以说可以用背包DP搞一搞。
代码：
/**************************************************************
Problem: 1004
User: danihao123
Language: C++
Result: Accepted
Time:156 ms
Memory:3172 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
#include <vector>
int sr, sb, sg, m, p;
int pow_mod(int a, int b) {
int ans = 1, res = a;
while(b) {
if(1 & b) ans = (ans * res) % p;
res = (res * res) % p; b >>= 1;
}
return ans;
}
int inv(int x) {
return pow_mod(x % p, p - 2);
}
int d[65][21][21][21];
std::vector<int> len;
int dp() {
int n = len.size();
d[0][0][0][0] = 1;
for(int i = 1; i <= n; i ++) {
int l = len[i - 1];
for(int j = 0; j <= sr; j ++) {
for(int k = 0; k <= sb; k ++) {
for(int t = 0; t <= sg; t ++) {
d[i][j][k][t] = 0;
if(j >= l) d[i][j][k][t] += d[i - 1][j - l][k][t];
if(k >= l) d[i][j][k][t] += d[i - 1][j][k - l][t];
if(t >= l) d[i][j][k][t] += d[i - 1][j][k][t - l];
d[i][j][k][t] %= p;
}
}
}
}
return d[n][sr][sb][sg];
}
int next[65]; bool vis[65];
int main() {
scanf("%d%d%d%d%d", &sr, &sb, &sg, &m, &p);
int n = sr + sb + sg;
int ans = 0;
for(int i = 1; i <= n; i ++) {
len.push_back(1);
}
ans = dp();
for(int i = 1; i <= m; i ++) {
memset(vis, 0, sizeof(vis));
for(int i = 1; i <= n; i ++) {
scanf("%d", &next[i]);
}
len.clear();
for(int i = 1; i <= n; i ++) {
if(!vis[i]) {
int p = i, cnt = 0;
do {
vis[p] = true; cnt ++;
p = next[p];
} while(p != i);
len.push_back(cnt);
}
}
ans = (ans + dp()) % p;
}
ans = (ans * inv(m + 1)) % p;
printf("%d\n", ans);
return 0;
}
[BZOJ 3640]JC的小苹果
逆矩阵文明，，，
很显然我们可以定义一个状态\(f[i][j]\)表示当前血量为\(i\)走到\(j\)的概率，然后肥肠爆歉这个东西没法DP（可能会有的点的伤害为0，这样可以凿出来环）。考虑高消，这个东西有个好处是不可能从血量低的到血量高的状态，所以可以从大到小枚举血量，这样各层事独立的，复杂度比直接高消降低了很少。可惜复杂度为\(O(sn^3)\)（设血量为\(s\)），会T掉。
考虑转移的过程，转移时等价于解这样一个方程：
\[
\begin{aligned}
a_{11}x_{1} + a_{12}x_{2} + \ldots + a_{1n}x_{n} &= c_1\\
a_{21}x_{1} + a_{22}x_{2} + \ldots + a_{2n}x_{n} &= c_2\\
&\ldots\\
a_{n1}x_{1} + a_{n2}x_{2} + \ldots + a_{nn}x_{n} &= c_n
\end{aligned}
\]
其中的未知数\(x\)事我们要求的东西，\(c\)表示从高血量状态转移过来的概率（这个可以视作常数）。根据友矩阵那一套理论，这一系列方程等价于以下等式：
\[
\begin{bmatrix}
a_{11} & a_{12} & \ldots & a_{1n}\\
a_{21} & a_{22} & \ldots & a_{2n}\\
\vdots & \vdots & \ddots & \vdots\\
a_{n1} & a_{n2} & \ldots & a_{nn}
\end{bmatrix}
\begin{bmatrix}
x_{1}\\ x_{2} \\ \vdots \\ x_{n}
\end{bmatrix}
=
\begin{bmatrix}
c_1\\ c_2\\ \vdots\\ c_n
\end{bmatrix}
\]
不妨将之简记为\(AB = C\)，其中我们只有\(B\)不知道，要求的也是\(B\)。我们除一下就可以得到\(B = A^{-1}C\)，然后我们还发现每一层的\(A\)都是一样的，所以每一层的\(A^{-1}\)也都是一样的，预处理即可。这样转移部分的复杂度就变成了\(O(sn^2)\)（矩阵乘法在这里事方阵乘列向量）。
至于逆矩阵的求法，我们知道对矩阵做初等变化也就等价于乘上另一个矩阵。因此，我们将一个矩阵\(A\)用类似于高消的手段消为单位阵\(I\)，所做的初等变换也就等价于乘上\(A^{-1}\)。我们对一个单位阵\(I\)作用上一样的操作，也就等于给这个单位阵乘上了\(A^{-1}\)，这样我们就得到了\(A^{-1}\)。
代码：
/**************************************************************
Problem: 3640
User: danihao123
Language: C++
Result: Accepted
Time:8708 ms
Memory:13416 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cassert>
#include <cmath>
#include <algorithm>
#include <functional>
#include <utility>
#include <vector>
#include <queue>
#include <set>
#include <map>
int n, m, hp;
const int maxn = 151, maxm = 5005;
const int maxh = 10005;
typedef double R;
typedef R Mat[maxn][maxn];
Mat D, F;
void calc_inv() {
for(int i = 1; i <= n; i ++) {
F[i][i] = 1;
}
for(int i = 1; i <= n; i ++) {
int r = i;
for(int j = i + 1; j <= n; j ++) {
if(fabs(D[j][i]) > fabs(D[r][i])) {
r = j;
}
}
// assert(r > -1);
if(r != i) {
for(int j = 1; j <= n; j ++) {
std::swap(D[r][j], D[i][j]);
std::swap(F[r][j], F[i][j]);
}
}
R bs = D[i][i];
for(int j = 1; j <= n; j ++) {
D[i][j] /= bs; F[i][j] /= bs;
}
for(int k = 1; k <= n; k ++) {
if(k != i) {
R rate = D[k][i];
for(int j = 1; j <= n; j ++) {
D[k][j] -= rate * D[i][j];
F[k][j] -= rate * F[i][j];
}
}
}
}
}
void matrix_mul(const Mat &A, const Mat &B, int a, int b, int c, Mat &res) {
static Mat C;
for(int i = 1; i <= a; i ++) {
for(int j = 1; j <= c; j ++) {
C[i][j] = 0;
}
}
for(int i = 1; i <= a; i ++) {
for(int j = 1; j <= c; j ++) {
for(int k = 1; k <= b; k ++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
for(int i = 1; i <= a; i ++) {
for(int j = 1; j <= c; j ++) {
res[i][j] = C[i][j];
}
}
}
int first[maxn], deg[maxn];
int next[maxm << 1], to[maxm << 1];
int gcnt = 0;
void add_edge(int u, int v) {
gcnt ++;
next[gcnt] = first[u]; first[u] = gcnt;
to[gcnt] = v;
}
void ins_edge(int u, int v) {
deg[u] ++;
add_edge(u, v);
if(u != v) {
deg[v] ++;
add_edge(v, u);
}
}
int atk[maxn];
R f[maxh][maxn];
R solve() {
for(int i = 1; i <= n; i ++) {
D[i][i] = 1.00;
}
for(int i = 1; i < n; i ++) {
for(int j = first[i]; j; j = next[j]) {
int v = to[j];
if(!atk[v]) {
D[v][i] -= 1.00 / (R(deg[i]));
}
}
}
calc_inv();
R ans = 0; static Mat C;
f[hp][1] = 1.00;
for(int i = hp; i >= 1; i --) {
for(int j = 1; j <= n; j ++) {
C[j][1] = f[i][j];
}
matrix_mul(F, C, n, n, 1, C);
for(int j = 1; j < n; j ++) {
for(int e = first[j]; e; e = next[e]) {
int v = to[e];
if(atk[v] > 0 && i - atk[v] > 0) {
f[i - atk[v]][v] += C[j][1] / (R(deg[j]));
}
}
}
ans += C[n][1];
}
return ans;
}
int main() {
scanf("%d%d%d", &n, &m, &hp);
for(int i = 1; i <= n; i ++) {
scanf("%d", &atk[i]);
}
for(int i = 1; i <= m; i ++) {
int u, v; scanf("%d%d", &u, &v);
ins_edge(u, v);
}
printf("%.8lf\n", solve());
return 0;
}
[BZOJ 5358][Lydsy1805月赛]口算训练
后几页有我会做的题……很好，我很安详，，，
考虑将所有数质因数分解。那么询问\([l, r]\)中所有数的积是否是\(d\)的倍数的本质就是对于\(d\)的每一类质因子，询问区间中该类质因子的指数之和是否不小于\(d\)中的。
考虑到数的范围都很小（不超过100000），我们可以先线性筛预处理，这样一次分解质因数的复杂度就降为了\(O(\log n)\)。至于维护区间每类质因子的指数和这种事……就用主席树处理吧。
代码：
/**************************************************************
Problem: 5358
User: danihao123
Language: C++
Result: Accepted
Time:2804 ms
Memory:68408 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
const int maxn = 100005;
int minp[maxn], mint[maxn], minh[maxn];
int prm[maxn], pcnt;
bool vis[maxn];
void sieve() {
const int N = 100000;
vis[1] = true; pcnt = 0;
for(int i = 2; i <= N; i ++) {
if(!vis[i]) {
prm[++ pcnt] = i;
minp[i] = pcnt; mint[i] = 1;
minh[i] = i;
}
for(int j = 1; j <= pcnt && i * prm[j] <= N; j ++) {
int v = i * prm[j];
vis[v] = true; minp[v] = j;
if(i % prm[j] == 0) {
mint[v] = mint[i] + 1; minh[v] = minh[i] * prm[j];
break;
} else {
mint[v] = 1; minh[v] = prm[j];
}
}
}
}
const int bufsiz = 64 * 1024 * 1024;
char buf[bufsiz]; char *cur;
void init_buf() {
cur = buf;
}
void *alloc(size_t size) {
if(cur + size - buf > bufsiz) {
return malloc(size);
} else {
char *ret = cur; cur += size;
return ret;
}
}
struct Node {
int v; Node *lc, *rc;
};
Node *build_tree(int L, int R) {
Node *ret = (Node*)alloc(sizeof(Node));
ret -> v = 0;
if(L == R) {
ret -> lc = ret -> rc = NULL;
} else {
int M = (L + R) / 2;
ret -> lc = build_tree(L, M);
ret -> rc = build_tree(M + 1, R);
}
return ret;
}
Node *add_p(Node *o, int L, int R, int p, int v) {
Node *ret = (Node*)alloc(sizeof(Node));
ret -> v = o -> v;
ret -> lc = o -> lc; ret -> rc = o -> rc;
if(L == R) {
ret -> v += v;
} else {
int M = (L + R) / 2;
if(p <= M) {
ret -> lc = add_p(ret -> lc, L, M, p, v);
} else {
ret -> rc = add_p(ret -> rc, M + 1, R, p, v);
}
}
return ret;
}
int query(Node *o, int L, int R, int p) {
if(L == R) {
return o -> v;
} else {
int M = (L + R) / 2;
if(p <= M) {
return query(o -> lc, L, M, p);
} else {
return query(o -> rc, M + 1, R, p);
}
}
}
Node *rt[maxn];
void solve() {
init_buf();
int n, m; scanf("%d%d", &n, &m);
rt[0] = build_tree(1, pcnt);
for(int i = 1; i <= n; i ++) {
rt[i] = rt[i - 1];
int x; scanf("%d", &x);
while(x > 1) {
int p = minp[x], t = mint[x];
rt[i] = add_p(rt[i], 1, pcnt, p, t);
x /= minh[x];
}
}
while(m --) {
int l, r, d; scanf("%d%d%d", &l, &r, &d);
bool ok = true;
while(d > 1) {
int p = minp[d], t = mint[d];
int st = query(rt[r], 1, pcnt, p) - query(rt[l - 1], 1, pcnt, p);
if(st < t) {
ok = false; break;
}
d /= minh[d];
}
if(ok) {
puts("Yes");
} else {
puts("No");
}
}
}
int main() {
sieve();
int T; scanf("%d", &T);
while(T --) {
solve();
}
return 0;
}
[BZOJ 4403]序列统计
老早做的题……
一看单调不降就想去用差分……但差分不好推（比下面的颓法要多一步……）。其实我们发现，只要给\([L, R]\)里每种整数分配出现次数，原序列就可以唯一确定了。
因此我们把\([L, R]\)中每个整数的出现次数当做一个变量，他们加起来应该等于一个\([1, n]\)中的整数。用隔板法很容易退出来式子是（令\(l = R - L + 1\)）：
\[\sum_{i = 1}^n \binom{i + l - 1}{l - 1}\]
看起来玩不动了……但是我们给式子加上一个\(\binom{l}{l}\)（其实就是1），然后我们会发现式子可以用杨辉三角的递推式合并成一个组合数，即\(\binom{n + l}{l}\)。然后求这个组合数还需要Lucas定理……
代码：
/**************************************************************
Problem: 4403
User: danihao123
Language: C++
Result: Accepted
Time:908 ms
Memory:16448 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <utility>
typedef long long ll;
const ll ha = 1000003LL;
const int maxn = 1000003;
ll pow_mod(ll a, ll b) {
ll ans = 1LL, res = a % ha;
while(b) {
if(1LL & b) ans = (ans * res) % ha;
res = (res * res) % ha;
b >>= 1;
}
return ans;
}
ll f[maxn], g[maxn];
void process() {
f[0] = 1LL;
for(int i = 1; i < maxn; i ++) {
f[i] = (f[i - 1] * (ll(i))) % ha;
}
g[maxn - 1] = pow_mod(f[maxn - 1], ha - 2LL);
for(int i = maxn - 2; i >= 0; i --) {
g[i] = (g[i + 1] * (ll(i + 1))) % ha;
}
}
ll C(int a, int b) {
if(a < b) return 0LL;
return (((f[a] * g[b]) % ha) * g[a - b]) % ha;
}
ll calc(int a, int b) {
if(a < b) return 0LL;
if(b == 0) return 1LL;
if(a < ha && b < ha) {
return C(a, b);
} else {
int as = a % ha, bs = b % ha;
return (C(as, bs) * calc(a / ha, b / ha)) % ha;
}
}
int main() {
process();
int T; scanf("%d", &T);
while(T --) {
int n, l, r; scanf("%d%d%d", &n, &l, &r);
int len = r - l + 1;
printf("%lld\n", (calc(n + len, len) - 1LL + ha) % ha);
}
return 0;
}
[BZOJ 3601]一个人的数论
本来想做JZPKIL的……但卡常太恶心了
上来先颓式子：
\[
\DeclareMathOperator{\id}{id}
\begin{aligned}
f_d(n) &= \sum_{i = 1}^n [(i, n) = 1]i^d\\
&= \sum_{i = 1}^n i^d\sum_{k | n, k | i}\mu(k)\\
&= \sum_{k | n}\mu(k)\sum_{i = 1}^{\frac{n}{k}} (ik)^d\\
&= \sum_{k | n}\mu(k)k^d h_d(\frac{n}{d})\\
&= ((\mu\cdot\id^k)\ast h_d)
\end{aligned}
\]
其中\(h_d\)表示指数为\(d\)时的等幂求和函数。然后我们可能会想，如果这个函数是积性函数，那么我们可以对质因数分解的每一种质因数都求一遍答案，而由于\(\mu\)的存在（对于质数的幂\(p^k\)的因数，只有\(1\)和\(p\)的莫比乌斯函数值不为0），所以需要处理的情况只有两种。
不过很可惜，\(h_d\)显然不是积性函数，所以最后的整个卷积也不是积性函数。但是我们注意到\(h_d\)事一个\(d + 1\)次多项式，所以我们可以把多项式每一项都当成单独的函数，然后每个单独函数卷出来都是一个积性函数。至于求多项式每一项的系数，用伯努利数即可。
代码：
/**************************************************************
Problem: 3601
User: danihao123
Language: C++
Result: Accepted
Time:220 ms
Memory:948 kb
****************************************************************/
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
typedef long long ll;
ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
ll ans = 1LL, res = a;
while(b) {
if(1LL & b) ans = (ans * res) % ha;
res = (res * res) % ha;
b >>= 1;
}
return ans;
}
ll inv(ll x) {
return pow_mod(x, ha - 2LL);
}
const int maxn = 105;
ll C[maxn][maxn];
void process_C() {
C[0][0] = 1LL;
for(int i = 1; i < maxn; i ++) {
C[i][0] = C[i][i] = 1LL;
for(int j = 1; j < i; j ++) {
C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % ha;
}
}
}
ll B[maxn];
void process_B() {
B[0] = 1LL;
for(int i = 1; i <= 101; i ++) {
B[i] = 1LL;
for(int j = 0; j < i; j ++) {
ll d = (((C[i][j] * B[j]) % ha) * inv(i - j + 1)) % ha;
B[i] = (B[i] - d + ha) % ha;
}
}
}
const int maxw = 1005;
typedef std::pair<ll, ll> pii;
pii P[maxw]; ll dv[maxw][3];
int d, w;
void process_dv() {
for(int i = 1; i <= w; i ++) {
ll p = P[i].first, t = P[i].second;
p %= ha;
dv[i][0] = pow_mod(p, t);
dv[i][1] = pow_mod(p, t - 1LL);
dv[i][2] = pow_mod(p, d);
}
}
ll calc(ll k, int z) {
ll ans = k;
for(int i = 1; i <= w; i ++) {
ll p = P[i].first, t = P[i].second;
p %= ha;
ll f1 = pow_mod(dv[i][0], z), f2 = pow_mod(dv[i][1], z);
ll v1 = (f1 - (dv[i][2] * f2) % ha + ha) % ha;
ans = (ans * v1) % ha;
}
#ifdef LOCAL
printf("Case (%lld, %d) : %lld\n", k, z, ans);
#endif
return ans;
}
int main() {
process_C(); process_B();
scanf("%d%d", &d, &w);
for(int i = 1; i <= w; i ++) {
scanf("%lld%lld", &P[i].first, &P[i].second);
}
process_dv();
ll ans = 0LL;
for(int i = 0; i <= d; i ++) {
ll g = calc((((C[d + 1][i] * B[i]) % ha) * inv(d + 1)) % ha, d + 1 - i);
ans = (ans + g) % ha;
}
printf("%lld\n", ans);
return 0;
}
[BZOJ 2321][BeiJing2011集训]星器
物理题可还行（其实我是想学势能方法，然后误入了……）
给每个星星（假设他在的位置是\((i, j)\)）定义势能\(E_p = i^2 + j^2\)，定义势函数\(\Phi(S)\)来表示状态\(S\)时所有星星势能的和。
然后我们发现，当两个星星互相吸引时（假设他们同行，坐标分别为\((i, j)\)和\((i, k)\)，且$j<k$），他们的坐标会变为\((i, j + 1)\)和\((i, k - 1)\)，势能的总减少量（也就是势函数的减小量）为\(2(k - j - 1)\)。
因此，整个过程中势函数的总减小量，就是答案的两倍。因此算出操作前后的势能做差即可……
代码：
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <utility>
typedef long long ll;
int main() {
int n, m; scanf("%d%d", &n, &m);
ll E = 0LL;
for(int i = 1; i <= n; i ++) {
for(int j = 1; j <= m; j ++) {
ll delta = i * i + j * j; ll x;
scanf("%lld", &x); delta *= x;
E += delta;
}
}
for(int i = 1; i <= n; i ++) {
for(int j = 1; j <= m; j ++) {
ll delta = i * i + j * j; ll x;
scanf("%lld", &x); delta *= x;
E -= delta;
}
}
E /= 2LL;
printf("%lld\n", E);
return 0;
}
