POJ - 2420 A Star not a Tree?（模拟退火/三分套三分）

传送门

---

题目大意

给出平面内的$n$个点找出费马点，即找到一个点使得该点到其他所有点的距离之和最小，输出这个最小距离。

解题思路

网上不少题解千篇一律的“假”模拟退火，比如$ACdreamer$这篇博客他自己在代码中就说这不算是模拟退火，因为当不采用“以一定概率接受比当前最优解更差的解”，就已经注定不是模拟退火，真正的模拟退火是个有概率依据的随机算法，具体参见下面代码：

//
// Created by Happig on 2020/11/9
//
#include <cstdio>
#include <math.h>
#include <algorithm>

using namespace std;

const double eps = 1e-14;
const double dinf = 1e300;
const double delta = 0.993;
const double T = 100;

struct Point {
    double x, y;

    Point(double a = 0, double b = 0) : x(a), y(b) {}
} p[105];

int n;
Point s;

int dcmp(double d) {
    if (fabs(d) < eps) return 0;
    return d > 0 ? 1 : -1;
}

double dis(Point a, Point b) {
    return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}

double f(Point cur) {
    double ans = 0;
    for (int i = 0; i < n; i++) {
        ans += dis(cur, p[i]);
    }
    return ans;
}

double SA() {
    double t = T, ans = dinf;
    Point tmp = s;
    while (t > eps) {
        Point nxt(tmp.x + (rand() * 2.0 - RAND_MAX) * t, tmp.y + (rand() * 2.0 - RAND_MAX) * t);
        double res = f(nxt);
        double del = res - ans;
        if (del < 0) {
            ans = res;
            s = tmp = nxt;
        } else if (exp((ans - res) / t) * RAND_MAX > rand()) {
            tmp = nxt;
        }
        t *= delta;
    }
    return ans;
}

void solve() {
    scanf("%d", &n);
    double sumx = 0, sumy = 0;
    for (int i = 0; i < n; i++) {
        scanf("%lf%lf", &p[i].x, &p[i].y);
        sumx += p[i].x, sumy += p[i].y;
    }
    s = Point(sumx / n, sumy / n);
    printf("%.0lf\n", SA());
}

int main() {
    //freopen("in.txt","r",stdin);
    //freopen("out.txt","w",stdout);
    //ios_base::sync_with_stdio(0),cin.tie(0),cout.tie(0);
    srand(100000007);
    srand(rand()), srand(rand());
    solve();
    return 0;
}

三分套三分

在洛谷的模拟退火推荐练习题中题解的一个思路启发了我。这个题可以三分来写。实际上对于费马点问题，有一个重要的结论（建议记住）：

当确定$x$坐标时，距离和函数关于$y$是一个单峰函数；确定$y$坐标时，距离和函数关于$x$是一个单峰函数。

于是我们可以先三分$x$再三分$y$，这样确定答案，之所以采用三分一百次的写法，是因为浮点误差不好控制，三分一百次的精度是最高的！

#include <cstdio>
#include <math.h>
#include <iostream>
#include <algorithm>

using namespace std;
const int maxn = 1010;
const double eps = 1e-10;

struct Point {
    double x, y;

    Point(double a = 0, double b = 0) : x(a), y(b) {}
} p[maxn];

int n;

double ly, ry;

double dis(Point a, Point b) {
    return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}

double f(double x, double y) {
    double ans = 0;
    Point cur(x, y);
    for (int i = 0; i < n; i++) {
        ans += dis(cur, p[i]);
    }
    return ans;
}

double tri_search(double x) {
    ly = 0, ry = 10000.0;
    double midl, midr;
    for (int i = 1; i <= 100; i++) {
        midl = ly + (ry - ly) / 3.0;
        midr = ry - (ry - ly) / 3.0;
        if (f(x, midl) > f(x, midr)) {
            ly = midl;
        } else ry = midr;
    }
    return f(x, midl);
}

int main() {
    scanf("%d", &n);
    for (int i = 0; i < n; i++) {
        scanf("%lf%lf", &p[i].x, &p[i].y);
    }
    double lx = 0, rx = 10000.0, midl, midr;
    for (int i = 1; i <= 100; i++) {
        midl = lx + (rx - lx) / 3.0;
        midr = rx - (rx - lx) / 3.0;
        if (tri_search(midl) > tri_search(midr)) {
            lx = midl;
        } else rx = midr;
    }
    printf("%.0lf\n", f(lx, ly));
    return 0;
}