HDU3694 四边形的费马点

Fermat Point in Quadrangle

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 65536/32768 K (Java/Others)
Total Submission(s): 1871    Accepted Submission(s): 325

Problem Description

In geometry the Fermat point of a triangle, also called Torricelli point, is a point such that the total distance from the three vertices of the triangle to the point is the minimum. It is so named because this problem is first raised by Fermat in a private letter. In the following picture, P ₀ is the Fermat point. You may have already known the property that:

Alice and Bob are learning geometry. Recently they are studying about the Fermat Point. 

Alice: I wonder whether there is a similar point for quadrangle.

Bob: I think there must exist one.

Alice: Then how to know where it is? How to prove?

Bob: I don’t know. Wait… the point may hold the similar property as the case in triangle. 

Alice: It sounds reasonable. Why not use our computer to solve the problem? Find the Fermat point, and then verify your assumption.

Bob: A good idea.

So they ask you, the best programmer, to solve it. Find the Fermat point for a quadrangle, i.e. find a point such that the total distance from the four vertices of the quadrangle to that point is the minimum.

 

Input

The input contains no more than 1000 test cases.

Each test case is a single line which contains eight float numbers, and it is formatted as below:

x ₁ y ₁ x ₂ y ₂ x ₃ y ₃ x ₄ y ₄

x _i, y _i are the x- and y-coordinates of the ith vertices of a quadrangle. They are float numbers and satisfy 0 ≤ x _i ≤ 1000 and 0 ≤ y _i ≤ 1000 (i = 1, …, 4).

The input is ended by eight -1.

 

Output

For each test case, find the Fermat point, and output the total distance from the four vertices to that point. The result should be rounded to four digits after the decimal point.

 

Sample Input

  
  

   
   0 0 1 1 1 0 0 1
1 1 1 1 1 1 1 1
-1 -1 -1 -1 -1 -1 -1 -1
  
  

 

Sample Output

  
  

   
   2.8284
0.0000
  
  

 

Source

2010 Asia Fuzhou Regional Contest

题意：给出四个点构成四边形（不按顺序给），求费马距离。
 
1.因为只有四个点，当这四个点可以构成凸四边形时（如图）：对角线形成的交点是费马点。就是点1，根据三角形的两边大于第三边可证明。
 
2.当不能构成凸多边形时，就是那个凹点是费马点（如图）：同理可证。

#define DeBUG
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <cstdlib>
#include <algorithm>
#include <vector>
#include <stack>
#include <queue>
#include <string>
#include <set>
#include <sstream>
#include <map>
#include <list>
#include <bitset>
#include <complex>
using namespace std ;
#define zero {0}
#define INF 0x3f3f3f3f
#define EPS 1e-6
#define TRUE true
#define FALSE false
typedef long long LL;
const double PI = acos(-1.0);
//#pragma comment(linker, "/STACK:102400000,102400000")
inline int sgn(double x)
{
    return fabs(x) < EPS ? 0 : (x < 0 ? -1 : 1);
}
#define N 100005
template<class T> T sqr(T x)//求平方
{
    return x * x;
}
// Point class
struct Point;
typedef Point Vec;
struct Point
{
    double x, y;
    Point () {}
    Point(double a, double b)
    {
        x = a;
        y = b;
    }
};

Vec operator + (const Vec &a, const Vec &b) //点加法
{
    return Vec(a.x + b.x, a.y + b.y);
}
Vec operator - (const Vec &a, const Vec &b) //点减法
{
    return Vec(a.x - b.x, a.y - b.y);
}
Vec operator * (const Vec &a, const double &p) //点与常数相乘
{
    return Vec(a.x * p, a.y * p);
}
Vec operator / (const Vec &a, const double &p) //点除以常数
{
    return Vec(a.x / p, a.y / p);
}
bool operator < (const Vec &a, const Point &b) //平面直角坐标系中左下方的为小
{
    return a.x < b.x || (sgn(a.x - b.x) == 0 && a.y < b.y);
}
bool operator == (const Vec &a, const Point &b) //点相等判断
{
    return sgn(a.x - b.x) == 0 && sgn(a.y - b.y) == 0;
}
Point p[10];
inline double ptDis(Point a, Point b)//点间距
{
    return sqrt(sqr(a.x - b.x) + sqr(a.y - b.y));
}
inline double dotDet(Vec a, Vec b)//点乘
{
    return a.x * b.x + a.y * b.y;
}
inline double crossDet(Vec a, Vec b)//叉乘
{
    return a.x * b.y - a.y * b.x;
}
inline double crossDet(Point o, Point a, Point b)//向量叉乘
{
    return crossDet(a - o, b - o);
}
// Polygon class
struct Poly//平面类
{
    vector<Point> pt;//保存平面对应的端点
    Poly()
    {
        pt.clear();
    }
    ~Poly() {}
    Poly(vector<Point> pt0): pt(pt0) {} //使用vector进行初始化
    Point operator [](int x) const//重载[]返回对应端点
    {
        return pt[x];
    }
    int size()//返回对应平面点数
    {
        return pt.size();
    }
    int isConvex()//0 不是多边形，1是凸多边形，-1不是凸多边形
    {
        int l = pt.size();
        double flag = crossDet(pt[0] - pt[l - 1], pt[1] - pt[0]);
        flag *= crossDet(pt[l - 1] - pt[l - 2], pt[0] - pt[l - 1]);
        for (int i = 1; i < l - 1; i++)
        {
            flag *= crossDet(pt[i] - pt[i - 1], pt[i + 1] - pt[i]);
        }
        if (sgn(flag) < 0)
            return -1;
        else if (sgn(flag) == 0)
            return 0;
        else
            return 1;
    }
} ;

bool cmp(const Point &p1, const Point &p2)//360度范围逆时针排序
{
    complex<double> c1(p1.x, p1.y), c2(p2.x, p2.y);
    return arg(c1) < arg(c2);
}

///逆时针顺序把极角排序 以p[0]（左下角点）为旋转变量
Point origin;
bool cmp3(const Point &a, const Point &b)
{
    int ret = crossDet(a - origin, b - origin);
    if (ret > 0)
        return true;
    if (ret == 0 && ptDis(origin, a) < ptDis(origin, b)) //保留凸包上的点，共线取最近
        return true;
    return false;
}
//引用上面排序例子
void PtAngleSort(Point *sp, Point *ep)
{
    int n = ep - sp;
    int k = 0;
    for (int i = 1; i < n; i++)
    {
        if (sp[i].x < sp[k].x || (sp[i].x == sp[k].x && sp[i].y < sp[k].y))
            k = i;
    }
    origin = sp[k];
    sort(sp, ep, cmp3);
}
int init()
{
    int flag = 0;
    for (int i = 0; i < 4; i++)
    {
        scanf("%lf%lf", &p[i].x, &p[i].y);
        if (p[i].x > 0 || p[i].y > 0)
            flag = 1;
    }
    // sort(p, p + 4, cmp);
    PtAngleSort(p, p + 4);
    // for (int i = 0; i < 4; i++)
    // {
    //     printf("%lf %lf\n", p[i].x, p[i].y);
    // }
    return flag;
}

int main()
{
#ifdef DeBUGs
    freopen("C:\\Users\\Sky\\Desktop\\1.in", "r", stdin);
    freopen("C:\\Users\\Sky\\Desktop\\dabiao.txt", "w", stdout);
#endif
    while (init())
    {
        std::vector<Point> v;
        for (int i = 0; i < 4; i++)
        {
            v.push_back(p[i]);
        }
        Poly poly(v);
        double ans = 0.0;
        if (poly.isConvex() == 1)
        {
            ans += ptDis(p[0], p[2]);
            ans += ptDis(p[1], p[3]);
        }
        else
        {
            p[4] = p[0];
            p[5] = p[1];
            bool flag = false;
            for (int i = 0; i < 4; i++)
            {
                if (crossDet(p[i + 1] - p[i], p[i + 2] - p[i + 1]) < 0)
                {
                    flag = true;
                    ans += ptDis(p[i + 1], p[(i + 3) % 4]);
                    ans += ptDis(p[i], p[i + 1]);
                    ans += ptDis(p[i + 1], p[i + 2]);
                    break;
                }
            }
            if (!flag)//这里有点坑爹啊，wa到死才加上这句就好了
            {
                ans += ptDis(p[0], p[2]);
                ans += ptDis(p[1], p[3]);
            }
        }
        printf("%.4lf\n", ans);
    }
    return 0;
}