题目地址:http://www.spoj.com/problems/DIST/
题目大意:给n个点(xi,yi),求(x0,y0),使得最小化sigma(|xi-x0|^pi+|yi-y0|^pi)。
算法讨论:因为xi和yi互不影响,因此可以分开计算。对于xi,先把xi排序,相邻2个xi之间的函数可以把绝对值去掉,这样的函数是一个一元三次函数。求一元三次函数的最值可以用导数,导数的零点就是原函数的极值点,然后再统计边界值即可。
Code:
#include <cmath>
#include <cstdio>
#include <cstring>
#include <utility>
#include <algorithm>
#define N 500000
#define eps 1e-9
#define inf 1000000000
#define x2(x) ((x)*(x))
#define x3(x) ((x)*(x)*(x))
using namespace std;
long double x0,x1;
long long f[4],ff[3];
pair<int,int> a[N+10];
int n,x[N+10],y[N+10],p[N+10];
inline void read(int &x){
x=0;
bool f=0;
char ch=getchar();
while (!(ch>='0' && ch<='9' || ch=='-')) ch=getchar();
if (ch=='-') f=1,ch=getchar();
while (ch>='0' && ch<='9') x=x*10+ch-'0',ch=getchar();
if (f) x=-x;
}
inline void solve(long double &x0,long double &x1){
long long a=ff[2],b=ff[1],c=ff[0],d=x2(b)-4*a*c;
if (a) if (d>=0) x0=(long double)(-b+sqrt(d))/2/a,x1=(long double)(-b-sqrt(d))/2/a;else x0=x1=1e9;else x0=x1=b?-(long double)c/b:1e9;
}
inline long double F(long double x){
return x3(x)*f[3]+x2(x)*f[2]+x*f[1]+f[0];
}
inline long double work(int x[]){
for (int i=1;i<=n;++i) a[i]=make_pair(x[i],p[i]);
sort(a+1,a+n+1);
memset(f,0,sizeof(f));
for (int i=1;i<=n;++i)
if (!a[i].second) f[0]++;//(xi-x)^0=1
else if (a[i].second==1) f[0]+=a[i].first,f[1]--;//(xi-x)^1=xi-x
else if (a[i].second==2) f[0]+=x2(a[i].first),f[1]-=2*a[i].first,f[2]++;//(xi-x)^2=xi^2-2xi*x+x^2
else f[0]+=x3(a[i].first),f[1]-=3*x2(a[i].first),f[2]+=3*a[i].first,f[3]--;//(xi-x)^3=xi^3-3xi^2*x+3xi*x^2-x^3
ff[2]=3*f[3],ff[1]=2*f[2],ff[0]=f[1];
long double ans=1e20;
a[0].first=-inf,a[n+1].first=inf;
for (int i=0;i<=n;++i){
solve(x0,x1);
if (a[i].first-x0<eps && x0-a[i+1].first<eps) ans=min(ans,F(x0));
if (a[i].first-x1<eps && x1-a[i+1].first<eps) ans=min(ans,F(x1));
ans=min(ans,F(a[i].first));
if (a[i+1].second==1) f[0]-=2*a[i+1].first,f[1]+=2;
else if (a[i+1].second==3) f[0]-=2*x3(a[i+1].first),f[1]+=6*x2(a[i+1].first),f[2]-=6*a[i+1].first,f[3]+=2;
ff[2]=3*f[3],ff[1]=2*f[2],ff[0]=f[1];
}
return ans;
}
int main(){
freopen("city.in","r",stdin);
freopen("city.out","w",stdout);
read(n);
for (int i=1;i<=n;++i) read(x[i]),read(y[i]),read(p[i]);
printf("%.3lf\n",(double)(work(x)+work(y)));
return 0;
}
By Charlie
Aug 22,2014