[bzoj3680] 吊打XXX：模拟退火 or 模拟力学情景

最新推荐文章于 2021-12-06 14:05:30 发布

chrt

最新推荐文章于 2021-12-06 14:05:30 发布

阅读量367

点赞数

分类专栏：近似算法模拟文章标签：模拟退火物理近似算法模拟

本文链接：https://blog.csdn.net/ruoruo_cheng/article/details/53081494

版权

近似算法同时被 2 个专栏收录

2 篇文章 0 订阅

订阅专栏

模拟

2 篇文章 0 订阅

订阅专栏

题意：n个（n<=10000）质量分别为 $m_i$ 的质点的坐标分别为 $(x_i, y_i)$ （数据均为整数），每个重物连着一根绳子，所有绳子有一个公共的绳结，绳结只能在一个与地面平行的平面上移动，质点不会落地，求绳结最终停留的坐标，保留小数点后三位。 $-100000 \le x_i, y_i \le 100000$ 。

可以用模拟退火，更简单的方法是模拟：给绳结一个初始位置，求受力，作用一段时间，产生一段位移。
嗯，都叫模拟……看是想模拟热力学还是力学 QAQ

模拟退火：
http://blog.csdn.net/popoqqq/article/details/39322569
https://zhuanlan.zhihu.com/p/21277465

这个问题其实叫广义费马点：
http://www.matrix67.com/blog/archives/422

#include <cstdio>
#include <cmath>
#include <cstdlib>
using namespace std;
const int MAX_N = 10000, loop = 1000;
const double pre = 1e-3, init = 1e4, c = 0.993, pi = acos(-1), inf = 1e70;
int n;
struct Point {
    double x, y, w;
    Point(double x=0, double y=0, double w=0): x(x), y(y), w(w) {}
    Point operator+(Point b)
    {
        return Point(x+b.x, y+b.y);
    }
    friend Point operator*(double k, Point p)
    {
        return Point(k*p.x, k*p.y);
    }
} P[MAX_N], a;
double sum, e = inf;

inline double dist(const Point& a, const Point& b)
{
    return sqrt((a.x-b.x)*(a.x-b.x) + (a.y-b.y)*(a.y-b.y));
}

inline double J(const Point& p)
{
    double r = 0;
    for (int i = 0; i < n; ++i)
        r += dist(p, P[i]) * P[i].w;
    if (r < e) {
        a = p;
        e = r;
    }
    return r;
}

inline double r()
{
    return (double)rand()/RAND_MAX;
}

inline void move(Point& p, Point& q, double d)
{
    double theta = 2*pi*r();
    q = p + d*Point(cos(theta), sin(theta));
}

void SA()
{
    Point p, q;
    for (int i = 0; i < n; ++i)
        p = p + P[i].w/sum*P[i];
    double last = J(p);
    for (double T = init; T >= pre; T *= c) {
        move(p, q, T);
        double now = J(q), dE = last - J(q);
        if (dE > 0 || r() < exp(dE/T)) {
            p = q;
            last = now;
        }
    }
    for (int i = 0; i < loop; ++i) {
        move(a, q, pre);
        J(q);
    }
}

int main()
{
    srand(2333333);
    scanf("%d", &n);
    for (int i = 0; i < n; ++i) {
        scanf("%lf %lf %lf", &P[i].x, &P[i].y, &P[i].w);
        sum += P[i].w;
    }
    SA();
    printf("%.3f %.3f\n", a.x, a.y);
    return 0;
}