题解 NKOJ2929 【[THUSC2014] 函数求解】

\(f(x) = x\) 显然是一组合法的解，但显然不一定优。

考虑利用特殊值给出更多性质：

• \(x = y\) 等价于 \(f(x) = f(y)\)。

证明：当 \(f(x) = f(y)\)，有 \(f(f(x) \times 1^m) = f(f(y) \times 1^m)\)，进而 \(x f(1)^m = y f(1)^m\)，可得 \(x = y\)。

由此我们可得 \(f\) 是一个双射。

接下来我们尝试探求 \(f(1)\) 的值。一个 naive 的想法是令 \(f(1) = 1\)，但这真的成立吗？

事实上这在 \(m \neq 2\) 时是成立的。证明：

• 考虑 \(f\) 的多层嵌套，如 \(f(f(f(f(1))))\)。
• 考虑一层一层拆括号，有 \(f(f(f(f(1)))) = f(f(f(f(1))) \times 1^m) = f(f(1)) \times f(1)^m = f(1)^{2m}\)。
• 考虑先拆里面一层，有 \(f(f(f(f(1)))) = f(f(f(f(1)) \times 1^m)) = f(f(1) \times f(1)^m) = f(f(1))^m = f(1)^{m^2}\)。
• 于是有 \(f(1) = 1\) 或 \(2m = m^2\)，则当 \(m \neq 2\) 时 \(f(1) = 1\)。

先讨论 \(m \neq 2\) 的情况，此时我们有如下结论：

• \(f(f(x)) = x\)。

证明：令 \(t = 1\) 即证。

• \(f(s \times t^m) = f(s) f(t)^m\)。

证明：令 \(s \leftarrow f(s)\) 即证。

• \(f(x^m) = f(x)^m\)。

证明：令 \(s = 1\) 即证。

• \(f\) 为完全积性函数。

证明：\(f(xy)^m = f(1 \times (xy)^m) = f((xy)^m) = f(x^m) f(y)^m = (f(x)f(y))^m\)。

于是我们的问题转化为：给每个质数 \(p\) 匹配一个质数 \(f(p)\)，答案为 \(\displaystyle\sum_{i = 1}^d q_i \ln f(p_i)\) 的最小值。

此时我们有两种情况：自己匹配自己或互相匹配。考虑从两者分别自己匹配自己到互相匹配的变化量，则问题转化为一般图最小权匹配。

注意到我们所用到的质数只有 \(\{p_i\}\) 和除此之外的前 \(d\) 小质数，则直接跑带权带花树的时间复杂度为 \(O(d^4)\)。想必能过！

但事实上直接把它当成二分图置于两侧跑费用流就是对的。

注意到我们要求的是最小费用流，则费用 \(\geq 0\) 的边不会被加上。

考虑存在的三条无向边 \((i, j), (i, k), (j, k)\)，我们已知 \(w(x, y) = (\ln_{p_x} - \ln_{p_y})(q_y - q_x)\)，则可知 \(w(i, j) > w(i, k) + w(j, k)\)，于是我们一定不会让一个点被重复选进匹配。

当 \(m = 2\)，\(f(1)\) 并不一定是 \(1\)，但感性理解一下 \(f(1) = 1\) 确实很优即可。

时间复杂度为 \(O(\operatorname{MF}(d, d^2))\)。记得在跑 SPFA 的时候判 eps。

代码：

#include <iostream>
#include <queue>
#include <cstdio>
#include <cmath>

using namespace std;

typedef struct {
	int nxt;
	int end;
	int dis;
	double cost;
} Edge;

const int N = 2e3, M = 400 + 7, K = 80800 + 7;
const double eps = 1e-7;
int cnt = 1;
int prime[N + 7], p[M], q[M], head[M], pre_dot[M], pre_edge[M];
double ln[N + 7], dis[M];
bool isnp[N + 7], mark[N + 7], vis[M];
Edge edge[K];
queue<int> que;

inline void init1(){
	int cnt = 0;
	for (register int i = 1; i <= N; i++){
		ln[i] = log(i);
	}
	isnp[0] = isnp[1] = true;
	for (register int i = 2; i <= N; i++){
		if (!isnp[i]) prime[++cnt] = i;
		for (register int j = 1; j <= cnt && i * prime[j] <= N; j++){
			isnp[i * prime[j]] = true;
			if (i % prime[j] == 0) break;
		}
	}
}

inline void init2(int n){
	for (register int i = 0; i <= n; i++){
		dis[i] = 1e9;
	}
}

inline void add_edge(int start, int end, int dis, double cost){
	cnt++;
	edge[cnt].nxt = head[start];
	head[start] = cnt;
	edge[cnt].end = end;
	edge[cnt].dis = dis;
	edge[cnt].cost = cost;
}

inline void spfa(int start){
	dis[start] = 0.0;
	vis[start] = true;
	que.push(start);
	while (!que.empty()){
		int cur = que.front();
		que.pop();
		vis[cur] = false;
		for (register int i = head[cur]; i != 0; i = edge[i].nxt){
			if (edge[i].dis != 0){
				int x = edge[i].end;
				double y = dis[cur] + edge[i].cost;
				if (dis[x] - y > eps){
					dis[x] = y;
					pre_dot[x] = cur;
					pre_edge[x] = i;
					if (!vis[x]){
						vis[x] = true;
						que.push(x);
					}
				}
			}
		}
	}
}

inline double mincost(int n, int start, int end){
	double ans = 0.0;
	while (true){
		init2(n);
		spfa(start);
		if (dis[end] >= 0.0) break;
		int flow = 0x7fffffff;
		for (register int i = end; i != start; i = pre_dot[i]){
			flow = min(flow, edge[pre_edge[i]].dis);
		}
		for (register int i = end; i != start; i = pre_dot[i]){
			edge[pre_edge[i]].dis -= flow;
			edge[pre_edge[i] ^ 1].dis += flow;
		}
		ans += flow * dis[end];
	}
	return ans;
}

int main(){
	int m, d, end;
	double base = 0.0;
	cin >> m >> d;
	init1();
	for (register int i = 1; i <= d; i++){
		cin >> p[i] >> q[i];
		mark[p[i]] = true;
		base += ln[p[i]] * q[i];
	}
	for (register int i = 1, j = 1; i <= d; i++, j++){
		while (mark[prime[j]]) j++;
		p[i + d] = prime[j];
	}
	d *= 2;
	end = d * 2 + 1;
	for (register int i = 1; i <= d; i++){
		int i_ = i + d;
		add_edge(0, i, 1, 0.0);
		add_edge(i, 0, 0, 0.0);
		for (register int j = 1; j <= d; j++){
			double delta = (ln[p[i]] - ln[p[j]]) * (q[j] - q[i]);
			if (delta < 0.0){
				int j_ = j + d;
				add_edge(i, j_, 1, delta);
				add_edge(j_, i, 0, -delta);
			}
		}
		add_edge(i_, end, 1, 0.0);
		add_edge(end, i_, 0, 0.0);
	}
	printf("%.4lf", base + mincost(end, 0, end) / 2.0);
	return 0;
}