模拟退火求二维费马点

一.概念引入

        AC过后，突然想起费马点和最小覆盖圆的圆心有关系吗？答案是没关系，如下图：费马点肯定在等腰梯形内，不是在圆心。

        模拟退火算法与初始值无关，算法求得的解与初始解状态S(是算法迭代的起点)无关；模拟退火算法具有渐近收敛性，已在理论上被证明是一种以概率收敛于全局最优解的全局优化算法；模拟退火算法具有并行性。

        三角形费马点:如果三角形ABC中有一个角的角度大于等于120, 那么该角的顶点就是费马点；如果三角形ABC的三个角的角度都小于120, 那么该三角形的费马点P为满足条件angle APB=angle BPC=angle CPA ，所以三角形的费马点也称为三角形的等角中心。

        牛顿迭代：令f(x,y)=Sigma( sqrt((x-xi)^2 + (y-yi)^2) ),本题要求f(x,y)的最小值，此函数是凸的。可验证极小值唯一，且当p0,p1..pn-1不在一条线上时，极小值点唯一，于是只需求得f'x=f'y=0的点即可。用牛顿迭代法解此二元方程组。

        塞瓦定理:设O是△ABC内任意一点，AO、BO、CO分别交对边于D、E、F，则 BD/DC*CE/EA*AF/FB=1。

        三角形内角平分线定理：今天给大二孩子将计算几何求内切圆圆心用到了。

二.算法实现

        以POJ2420为例，直接去AC吧。

        模拟退火（Simulated Annealing），题目：就是求多边形的费马点，输出最小的距离。讨论群里有人说应该说牛顿迭代是……所有大学生应该想到的……枚举是……高中生想到的……顿时无语。

import java.util.Scanner;

public class POJ2420 {

  static Point[] p;

  static int n;

  public static void main(String[] args) {

    double x;

    double y;

    Scanner sc = new Scanner(System.in);

    while(sc.hasNext()) {

      n = sc.nextInt();

      p = new Point[n];

      for(int i=0; i<n; i++) {

        p[i] = new Point();

      }

      for(int i=0; i<n; i++) {

        x = sc.nextDouble();

        y = sc.nextDouble();

        p[i] = new Point(x,y);

      }

      double ans = solve();

      System.out.println((int)ans);

    }

  }

  private static double solve() {

    Point pre = new Point();

    Point cur = new Point();

    /*

     * 终于找到错了，发现"cur = pre"的话，cur和pre完全一样，

     * 调试里id都是34，内存一样，

     * 我把其他的赋值全部搞成属性赋值就AC了

     * 坑了一天

     */

    double dis = getAllDistance(cur);

    double step = 10000;

    double r = 0.5;

    double[][] d = {{0,1},{0,-1},{-1,0},{1,0}};

    Point node = new Point();

    while(Double.compare(step, 0.2)>0) {

      boolean ok = true;

      while(ok) {

        ok = false;

        /*

         * 下面的赋值是必须的

         * 因为如果下面的for循环没有更新node，

         * 那么for循环后的赋值就没法进行

         */

        node.x = cur.x;

        node.y = cur.y;

        for(int i=0; i<4; i++) {

          //四次循环中cur不能 变

          pre.x = d[i][0]*step + cur.x;

          pre.y = d[i][1]*step + cur.y;

          //由于坐标大于0，所以小于0的坐标定然不是最小值

          if(pre.x<0||pre.y<0) {

            continue;

          }

          double temp = getAllDistance(pre);

//          System.out.println("step:"+step +"  " + 

//          "dis:"+dis +"   "+"temp:"+temp+

//          " "+cur.x+"  "+cur.y+"---"+pre.x+"  "+pre.y);

          if(dis>temp) {

            dis = temp;

            node.x = pre.x;

            node.y = pre.y;

            ok = true;

          }

        }

        cur.x = node.x;

        cur.y = node.y;

      }

      step *= r;

    }

    return Math.floor(dis+0.5);

    //return dis;

  }

  private static double getAllDistance(Point cur) {

    double sum = 0;

    for(int i=0; i<n; i++) {

      sum += getDistance(cur,p[i]);

    }

    return sum;

  }

  private static double getDistance(Point a, Point b) {

    return Math.hypot(a.x-b.x, a.y-b.y);

  }

}

class Point {

  double x;

  double y;

  public Point() {

    this.x = 0;

    this.y = 0;

  }

  public Point(double x, double y) {

    super();

    this.x = x;

    this.y = y;

  }

}