hdu 4717 The Moving Points(三分法)

// Time 78 ms; Memory 1316K  

#include<iostream>  
#include<cstdio>  
#include<cmath>  
#define maxn 50000  
#define maxm 305  
#define sqr(x) ((x)*(x))  
#define eps 1e-5  
  
using namespace std;  
  
double a[maxn],b[maxn],c[maxn];  
int cas,cnt;  
  
int sig(double x)  
{  
    return (x>eps)-(x<-eps);  
}  
  
typedef struct point  
{  
    double x,y;  
    point(double xx=0,double yy=0):x(xx),y(yy){}  
}vector;  
  
point pt[maxm];  
vector vt[maxm];  
  
vector operator - (point p,point q)  
{  
    return vector(p.x-q.x,p.y-q.y);  
}  
double dot(vector p,vector q)  
{  
    return p.x*q.x+p.y*q.y;  
}  
double cal(double t)  
{  
    double tmp,res=a[0]*t*t+b[0]*t+c[0];  
    for(int i=1;i<cnt;i++)  
    {  
        tmp=a[i]*t*t+b[i]*t+c[i];  
        if(res<tmp) res=tmp;  
    }  
    return res;  
}  
  
int main()  
{  
    int i,j,t,n;  
    double l,r,m1,m2,tmp;  
    vector v,w;  
    scanf("%d",&t);  
    while(t--)  
    {  
        scanf("%d",&n);  
        cnt=0;  
        l=r=0.0;  
        for(i=0;i<n;i++)  
            scanf("%lf%lf%lf%lf",&pt[i].x,&pt[i].y,&vt[i].x,&vt[i].y);  
        for(i=0;i<n;i++) for(j=i+1;j<n;j++)  
        {  
            v=vt[j]-vt[i],  
            w=pt[j]-pt[i],  
            a[cnt]=dot(v,v),  
            b[cnt]=2*dot(v,w),  
            c[cnt]=dot(w,w);  
            tmp=-b[cnt]/(2*a[cnt]);  
            if(r<tmp) r=tmp;  
            cnt++;  
        }  
        while(sig(r-l)>0)  
        {   
            m1=l+(r-l)/3;m2=r-(r-l)/3;  
            if(cal(m1)<cal(m2)) r=m2;  
            else l=m1;  
        }  
        printf("Case #%d: %.2lf %.2lf\n",++cas,l,sqrt(cal(l)));  
    }  
    return 0;  
}