HDU 2440、HDU 3694多边形费马点

费马点是指求到N个点距离和最小的点。通常要求的是简单多边形的费马点。

1、http://acm.hdu.edu.cn/showproblem.php?pid=2440   按照题意知道是一个简单的多边形即凸包，但给出的点并没有按照顺序的，所以需要自己先求出凸包，然后在用随机淬火求费马点。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<cmath>
#include<cstdlib>
#include<queue>
#include<stack>
#include<map>
#include<vector>
#include<algorithm>
#include<ctime>
using namespace std;
#define eps 1e-10

int Fabs(double d)
{
	if(fabs(d)<eps) return 0;
	else return d>0?1:-1;
}

struct point 
{
	double x,y;
}p[105],sta[105];
int oper[8][2]={0,1,0,-1,-1,0,1,0,1,1,1,-1,-1,1,-1,-1},top;

double x_multi(point p1,point p2,point p3)    
{  
    return (p2.x-p1.x)*(p3.y-p1.y)-(p3.x-p1.x)*(p2.y-p1.y);  
} 

double Dis(point p1,point p2)
{
	return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
}

bool cmp(point a,point b)
{
	if(Fabs(x_multi(p[0],a,b))>0) return 1;
	if(Fabs(x_multi(p[0],a,b))<0) return 0;
	if(Fabs(Dis(p[0],a)-Dis(p[0],b))<0)
		return 1;
	return 0;
}

void Graham(int n)
{
	int i,k=0,tot;
	for(i=1;i<n;i++)
		if((p[i].y<p[k].y)||((p[i].y==p[k].y)&&(p[i].x<p[k].x)))
			k=i;
	swap(p[0],p[k]);
	sort(p+1,p+n,cmp);
	
	tot=1;
    for(i=2;i<n;i++)
		if(Fabs(x_multi(p[i],p[i-1],p[0])))
			p[tot++]=p[i-1];
	p[tot++]=p[n-1];
	
	sta[0]=p[0],sta[1]=p[1];
	i=top=1;
	for(i=2;i<tot;i++)
	{ 
        while(top>=1&&Fabs(x_multi(p[i],sta[top],sta[top-1]))>=0)
		{ 
            if(top==0) break; 
            top--; 
        } 
        sta[++top]=p[i]; 
    }   
}

double allDis(int n,point f)
{
	double sum=0.0;
	int i;
	for(i=0;i<n;i++)
		sum+=Dis(sta[i],f);
	return sum;
}

point fermat(int n) //求费马点
{
	double step=0;
	int i,j;
	for(i=0;i<n;i++)
		step+=fabs(sta[i].x)+fabs(sta[i].y);
	point f;
	f.x=0,f.y=0;
	for(i=0;i<n;i++)
		f.x+=sta[i].x,f.y+=sta[i].y;
	f.x/=n,f.y/=n;
	point t;
	while(step>eps)
	{
		for(i=0;i<8;i++)
		{
			t.x=f.x+oper[i][0]*step;
			t.y=f.y+oper[i][1]*step;
			if(allDis(n,t)<allDis(n,f))
				f=t;
		}
		step/=2;
	}
	return f;
}

int main()
{
	int t,n,i;
	scanf("%d",&t);
	while(t--)
	{
		scanf("%d",&n);
		for(i=0;i<n;i++)
			scanf("%lf%lf",&p[i].x,&p[i].y);
		Graham(n);
		point ans=fermat(top+1);
		printf("%.0lf\n",allDis(top+1,ans));
		if(t>0) puts("");
	}
	return 0;
}

2、http://acm.hdu.edu.cn/showproblem.php?pid=3694 //求一个四边形的费马点，wrong了n次网上到处查才知道此题非常严谨，卡随机淬火算法。并且给出的四边形并不一定是凸四边形，所以需要讨论，如果是凸四边形，按照四边形的特性费马点就是对角线的交点，如果是凹的就是其中某一个顶点。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<cmath>
#include<cstdlib>
#include<queue>
#include<stack>
#include<map>
#include<vector>
#include<algorithm>
#include<ctime>
using namespace std;
#define eps 1e-8

int Fabs(double d)
{
    if(fabs(d)<eps) return 0;
    else return d>0?1:-1;
}

struct point 
{
    double x,y;
}p[10],sta[10];
int oper[8][2]={0,1,0,-1,-1,0,1,0,1,1,1,-1,-1,1,-1,-1},top;

double x_multi(point p1,point p2,point p3)    
{  
    return (p2.x-p1.x)*(p3.y-p1.y)-(p3.x-p1.x)*(p2.y-p1.y);  
} 

double Dis(point p1,point p2)
{
    return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
}

bool cmp(point a,point b)
{
    if(Fabs(x_multi(p[0],a,b))>0) return 1;
    if(Fabs(x_multi(p[0],a,b))<0) return 0;
    if(Fabs(Dis(p[0],a)-Dis(p[0],b))<0)
        return 1;
    return 0;
}

void Graham(int n)
{
    int i,k=0,tot;
    for(i=1;i<n;i++)
        if((p[i].y<p[k].y)||((p[i].y==p[k].y)&&(p[i].x<p[k].x)))
            k=i;
    swap(p[0],p[k]);
    sort(p+1,p+n,cmp);
    
    /*tot=1;//下面直接用顶点个数判断是否为凸包，所以这里不去共线点
    for(i=2;i<n;i++)
		if(Fabs(x_multi(p[i],p[i-1],p[0])))
			p[tot++]=p[i-1];
    p[tot++]=p[n-1];*/
    
    sta[0]=p[0],sta[1]=p[1];
    i=top=1;
    for(i=2;i<n;i++)
    { 
        while(top>=1&&Fabs(x_multi(p[i],sta[top],sta[top-1]))>=0)
        { 
            if(top==0) break; 
            top--; 
        } 
        sta[++top]=p[i]; 
    }   
}

/*double allDis(int n,point f)
{
    double sum=0.0;
    int i;
    for(i=0;i<n;i++)
        sum+=Dis(p[i],f);
    return sum;
}

point fermat(int n)
{
    double step=0;
    int i,j;
    for(i=0;i<n;i++)
        step+=fabs(sta[i].x)+fabs(sta[i].y);
    point f;
    f.x=0,f.y=0;
    for(i=0;i<n;i++)
        f.x+=sta[i].x,f.y+=sta[i].y;
    f.x/=n,f.y/=n;
    point t;
    while(step>1e-10)
    {
        for(i=0;i<8;i++)
        {
            t.x=f.x+oper[i][0]*step;
            t.y=f.y+oper[i][1]*step;
            if(allDis(n,t)<allDis(n,f))
                f=t;
        }
        step/=2;
    }
    return f;
}
*/
int main()
{
	int i,j;
    while(~scanf("%lf%lf%lf%lf%lf%lf%lf%lf",&p[0].x,&p[0].y,&p[1].x,&p[1].y,&p[2].x,&p[2].y,&p[3].x,&p[3].y))
    {
        if(p[0].x==-1&&p[0].y==-1&&p[1].x==-1&&p[1].y==-1&&p[2].x==-1&&p[2].y==-1&&p[3].x==-1&&p[3].y==-1)
            break;
        Graham(4);
        double ans;
        if(top==3)
			ans=Dis(sta[0],sta[2])+Dis(sta[1],sta[3]);//凸四边形就直接取对角线交点
        else 
        {
        	ans=1e50;
			double sum=0;
        	for(i=0;i<4;i++)
        	{
            	sum=0.0;
            	for(j=0;j<4;j++)
                	if(i!=j)
                    	sum+=Dis(p[i],p[j]);
             	ans=min(sum,ans);
        	}
		}
		printf("%.4lf\n",ans);
	}
    return 0;
}