FZU 1576【计算几何/费马点】

Oaiei居住在A城市，并且是这个城市建设的总设计师。最近有个问题一直困恼着他。A城市里有三个大型工厂，每个大型工厂每天都需要消耗大量的石油，现在城市里要建设一个石油中转站，从石油中转站到三个大型工厂都需要铺设石油管道。现在你的问题来了，应该如何建设这个石油中转站，使得石油中转站到三个大型工厂所需要铺设的石油管道线路最短，你能够帮助他吗？

设石油中转站B的坐标为(X,Y)，三个大型工厂的坐标分别为C(X1,Y1)，D(X2,Y2)，E(X3,Y3)，所输要铺设的石油管道线路总长为distance(B,C)+distance(B,D)+distance(B,E)，其中distance(A,B)表示A,B两点之间的欧几里德距离。

Input

第一行为一个整数C(1<=C<=200)，表示测试数据的组数。
以下C行，每行6个整数，分别表示C、D、E点 (X,Y)的坐标。注意：B点的选址可以与C、D、E点相重合。

Output

对于每个输入数据，输出一行两个数M1,M2，中间用一个空格分隔开。M1，M2分别四舍五入到小数点后两位，表示石油中转站的坐标。

Sample Input

1
0.00 0.00 1.00 0.00 0.00 1.00

Sample Output

0.21 0.21 

【模板】：http://blog.csdn.net/dingyaguang117/article/details/7216479

1.在一个三角形中，到3个顶点距离之和最小的点叫做这个三角形的费马点。

2.费马点计算方法:

(1)若三角形ABC的3个内角均小于120°，那么3条距离连线正好平分费马点所在的周角。所以三角形的费马点也称为三角形的等角中心。

(2)若三角形有一内角不小于120度，则此钝角的顶点就是距离和最小的点。 

3.如何计算等角中心呢？

做任意一条边的外接等边三角形，得到另一点，将此点与此边在三角形中对应的点相连

如此再取另一边作同样的连线，相交点即费马点

 

证明画几条辅助线就出来了~这里就不证明了

【代码】：

#include <iostream>
#include <math.h>
#include <windows.h>
#include <iomanip>
using namespace std;

struct Vec
{
    double x,y;
    Vec(double xx=0,double yy=0)
    {
        x=xx;
        y=yy;
    }
};

struct Point
{
    double x,y;
    Point(double xx=0,double yy=0)
    {
        x=xx;
        y=yy;
    }
};

double ddot(Vec A,Vec B)
{
    return A.x*B.x+A.y*B.y;
}
double getlen(Vec A)
{
    return sqrt(A.x*A.x+A.y*A.y);
}

double getlen(Point A,Point B)
{
    return sqrt((A.x-B.x)*(A.x-B.x)+(A.y-B.y)*(A.y-B.y));
}

bool moreThan120(double xa,double ya,double xb,double yb,double xc,double yc)
{
    Vec ab(xb-xa,yb-ya),ac(xc-xa,yc-ya);
    if (ddot(ab,ac)/getlen(ab)/getlen(ac) < -0.5)
    {
        return true;
    }
    return false;
}

inline void swap(double &a,double &b)
{
    double t;
    t=a;
    a=b;
    b=t;
}


Point getAnotherPoint(Point A,Point B,Point C)
{
    Point r,C1,C2;
    Vec AB(B.x-A.x,B.y-A.y);
    double len,len2;
    double sqrt3=sqrt(3.0);
    Vec AB2,crossAB,crossAB2;

    AB2.x=AB.x/2; AB2.y=AB.y/2;
    crossAB.x=AB2.y; crossAB.y=-AB2.x;
    crossAB2.x=-AB2.y; crossAB2.y=AB2.x;

    len=getlen(AB2);

    crossAB.x*=sqrt3; crossAB.y*=sqrt3;
    crossAB2.x*=sqrt3; crossAB2.y*=sqrt3;

    C1.x=A.x+AB2.x+crossAB.x;C1.y=A.y+AB2.y+crossAB.y;
    C2.x=A.x+AB2.x+crossAB2.x;C2.y=A.y+AB2.y+crossAB2.y;

    if (getlen(C,C1)<getlen(C,C2))
    {
        return C2;
    }else
        return C1;

}


/*
在平面内的两条直线AB,CD，求交点最直接的方法就是解下列的二元二次方程组：
Ax + (Bx - Ax)i = Cx + (Dx - Cx) j
Ay + (By - Ay)i = Cy + (Dy - Cy) j
交点是:
(Ax + (Bx - Ax) i, Ay + (By - Ay) i)
即：
  Ax + (AAx)i = Bx + (BBx) j
  Ay + (AAy)i = By + (BBy) j
  交点是:
(Ax + (AAx)i, Ay + (AAy)i)
*/
Point getCrossPoint(Point A,Point A1,Point B,Point B1)
{
    Point r;
    Vec AA(A1.x-A.x,A1.y-A.y),BB(B1.x-B.x,B1.y-B.y);
    double i,j,tmp,tmp2;
    double Ax=A.x,Ay=A.y,AAx=AA.x,AAy=AA.y,Bx=B.x,By=B.y,BBx=BB.x,BBy=BB.y;

    if (AAx==0)
    {
        j=(Ax-Bx)/BBx;
        i=(By+BBy*j-Ay)/AAx;
    }else if (BBx==0)
    {
        i=(Bx-Ax)/AAx;
    }else if (AAy==0)
    {
        j=(Ay-By)/BBy;
        i=(Bx-Ax-BBx*j)/AAx;
    }else if (BBy==0)
    {
        i=(By-Ay)/AAy;
    }
    else
    {
        tmp=AAx;
        tmp2=AAy;
        Ax*=AAy;AAx*=AAy;Bx*=AAy;BBx*=AAy;
        Ay*=tmp;AAy*=tmp;By*=tmp;BBy*=tmp;
        j=((Ax-Ay)-(Bx-By))/(BBx-BBy);
        i=(Bx+BBx*j-Ax)/AAx;
    }

    r.x=(Ax+AAx*i)/tmp2;
    r.y=(Ay+AAy*i)/tmp;
    return r;
}

int main()
{
    //freopen("cul.in9","r",stdin);
    int t;
    cin>>t;
    while(t--)
    {
        double xa,ya,xb,yb,xc,yc;
    Point C1,A1,R;


    cin>>xa>>ya>>xb>>yb>>xc>>yc;
    cout.setf(ios::fixed);

    if (moreThan120(xa,ya,xb,yb,xc,yc))
    {
        cout<<setprecision(2)<<xa<<" "<<ya<<endl;
    }else if (moreThan120(xb,yb,xa,ya,xc,yc))
    {
        cout<<setprecision(2)<<xb<<" "<<yb<<endl;
    }else if (moreThan120(xc,yc,xa,ya,xb,yb))
    {
        cout<<setprecision(2)<<xc<<" "<<yc<<endl;
    }else
    {
        C1=getAnotherPoint(Point(xa,ya),Point(xb,yb),Point(xc,yc));
        A1=getAnotherPoint(Point(xc,yc),Point(xb,yb),Point(xa,ya));

        R=getCrossPoint(Point(xa,ya),A1,Point(xc,yc),C1);
        cout<<setprecision(2)<<R.x<<" "<<R.y<<endl;
    }
    }
    //Sleep(1000000);
}