模拟退火学习笔记

此处我们以 [JSOI2004]平衡点 为例。
这题是叫你在一个二维平面上找出若干个力的平衡点。
由于不知道这个点在何处，所以考虑乱猜。
但猜完以后我们需要找更优的解，这是可以考虑随机一个。
这样的话，设猜的点坐标为 \(x,y\) ，随机的新点坐标为 \(x_1,y_1\) 。
则如果 \(x_1,y_1\) 比 \(x,y\) 更优的话，我们就接受这个新的点 \(x_1,y_1\) 为当前最优解。
如果 \(x_1,y_1\) 比不过 \(x,y\) ，为防止解就死在 \(x,y\) 周围的一个小范围内，我们还是要有一定概率接受 \(x_1,y_1\) 这个点的。但是随着时间的推移，我们越来越接近最优解时，接受一个不太优秀的点就要谨慎了。这就是所谓的“降温”操作。
等到温度降到 \(0\) 时，一切都结束了。

代码：

#include<stdio.h>
#include<time.h>
#include<stdlib.h>
#include<math.h>
#define reg register
#define ri reg int
#define rep(i, x, y) for(ri i = x; i <= y; ++i)
#define nrep(i, x, y) for(ri i = x; i >= y; --i)
#define DEBUG 1
#define delta 0.996
#define ll long long
#define il inline
#define max(i, j) (i) > (j) ? (i) : (j)
#define min(i, j) (i) < (j) ? (i) : (j)
#define read(i) io.READ(i)
#define print(i) io.WRITE(i)
#define push(i) io.PUSH(i)
struct IO {
#define MAXSIZE (1 << 20)
#define isdigit(x) (x >= '0' && x <= '9')
  char buf[MAXSIZE], *p1, *p2;
  char pbuf[MAXSIZE], *pp;
#if DEBUG
#else
  IO() : p1(buf), p2(buf), pp(pbuf) {}
  ~IO() {
    fwrite(pbuf, 1, pp - pbuf, stdout);
  }
#endif
  inline char gc() {
#if DEBUG
    return getchar();
#endif
    if(p1 == p2)
      p2 = (p1 = buf) + fread(buf, 1, MAXSIZE, stdin);
    return p1 == p2 ? ' ' : *p1++;
  }
  inline bool blank(char ch) {
    return ch == ' ' || ch == '\n' || ch == '\r' || ch == '\t';
  }
  template <class T>
  inline void READ(T &x) {
    register double tmp = 1;
    register bool sign = 0;
    x = 0;
    register char ch = gc();
    for(; !isdigit(ch); ch = gc())
      if(ch == '-') sign = 1;
    for(; isdigit(ch); ch = gc())
      x = x * 10 + (ch - '0');
    if(ch == '.')
      for(ch = gc(); isdigit(ch); ch = gc())
        tmp /= 10.0, x += tmp * (ch - '0');
    if(sign) x = -x;
  }
  inline void READ(char *s) {
    register char ch = gc();
    for(; blank(ch); ch = gc());
    for(; !blank(ch); ch = gc())
      *s++ = ch;
    *s = 0;
  }
  inline void READ(char &c) {
    for(c = gc(); blank(c); c = gc());
  }
  inline void PUSH(const char &c) {
#if DEBUG
    putchar(c);
#else
    if(pp - pbuf == MAXSIZE) {
      fwrite(pbuf, 1, MAXSIZE, stdout);
      pp = pbuf;
    }
    *pp++ = c;
#endif
  }
  template <class T>
  inline void WRITE(T x) {
    if(x < 0) {
      x = -x;
      PUSH('-');
    }
    static T sta[35];
    T top = 0;
    do {
      sta[top++] = x % 10;
      x /= 10;
    } while(x);
    while(top)
      PUSH(sta[--top] + '0');
  }
  template <class T>
  inline void WRITE(T x, char lastChar) {
    WRITE(x);
    PUSH(lastChar);
  }
} io;
int n;
struct node {
  int x, y, w;
} obj[2010];
double ansx, ansy, answ;
double calc(double x, double y) { /*x, y是X点的坐标*/
  double r = 0/*合力*/, dx/*横轴上的合力*/, dy/*纵轴上的合力*/;
  rep(i, 1, n) {
    dx = x - obj[i].x;
    dy = y - obj[i].y;
    r += sqrt(dx * dx + dy * dy) * obj[i].w; /*计算每一个质点产生的力*/
  }
  return r;
}
void anneal() { 
  double t = 3000; /*初始温度(可能要调)*/
  while(t > 1e-15) { 
    double rx = ansx + (rand() * 2 - RAND_MAX) * t;
    double ry = ansy + (rand() * 2 - RAND_MAX) * t;
    /*随机一个X点*/ 
    double rw = calc(rx, ry);
    /*计算这个点的答案*/ 
    if (rw - answ < 0) { /*当前答案更优,那么无条件接受*/
      ansx = rx;
      ansy = ry;
      answ = rw;
    } else if(exp((-rw + answ) / t) * RAND_MAX > rand()) { 
      ansx = rx;
      ansy = ry;
    }
    t *= delta; /*降温*/
  }
}
int main() {
  read(n);
  rep(i, 1, n) 
    read(obj[i].x), read(obj[i].y), read(obj[i].w), ansx += obj[i].x, ansy += obj[i].y;
  ansx /= n;
  ansy /= n;
  /*以平均数为初始答案*/
  answ = calc(ansx, ansy);
  while((double)clock() / CLOCKS_PER_SEC < 0.8) anneal(); /*模拟退火*/
  printf("%.3lf %.3lf\n", ansx, ansy);
  return 0;
}