D.Country Meow 最小球覆盖 三分套三分套三分 && 模拟退火

// 2019.10.3

// 练习题：2018 ICPC 南京现场赛

D Country Meow

题目大意

给定空间内 N 个点，求某个点到 N 个点的距离最大值的最小值。

 

思路

非常裸的最小球覆盖问题啊，即找到半径最小的球包含全部的点。

在最小圆覆盖问题上，可以使用随机增量法，这里没有四点确定球心的公式，所以板子失效了。

最小圆覆盖可以用三分套三分，这里空间有三维，假装证明得到在任意一维上都满足凸函数特性，那么再套一层维度三分就OK了。

 

AC代码

三分套三分套三分写法，复杂度O(n*log^3)

#include<iostream>

#include<cstdio>

#include<cmath>

#include<algorithm>

using namespace std;

const double eps = 1e-3;

struct Point {

	double x, y, z;

	Point() {

		x = y = z = 0.0;

	}

	Point(double xx, double yy, double zz) {

		x = xx, y = yy, z = zz;

	}

	Point operator-(const Point& p) {

		return Point(x-p.x, y-p.y, z-p.z);

	}

	double dis() {

		return sqrt(x*x+y*y+z*z);

	}

}pt[110];

int n;

double cal(double x, double y, double z) {

	double res = 0;

	for(int i=1;i<=n;i++) {

		res = max(res, (pt[i]-Point(x, y, z)).dis());

	}

	return res;

}

double cal2(double x, double y) {

	double res = 1e18;

	double l = -100000, r = 100000;

	while(r-l>eps) {

		double m1 = (r-l)/3 + l;

		double m2 = (r-l)/3*2 + l;

		double res1 = cal(x, y, m1), res2 = cal(x, y, m2);

		res = min(res, min(res1, res2));

		if(res1<res2) r = m2;

		else l = m1;

	}

	return res;

}

double cal3(double x) {

	double res = 1e18;

	double l = -100000, r = 100000;

	while(r-l>eps) {

		double m1 = (r-l)/3 + l;

		double m2 = (r-l)/3*2 + l;

		double res1 = cal2(x, m1), res2 = cal2(x, m2);

		res = min(res, min(res1, res2));

		if(res1<res2) r = m2;

		else l = m1;

	}

	return res;

}

int main() {

	cin>>n;

	for(int i=1;i<=n;i++) {

		scanf("%lf %lf %lf", &pt[i].x, &pt[i].y, &pt[i].z);

	}

	double res = 1e18;

	double l = -100000, r = 100000;

	while(r-l>eps) {

		double m1 = (r-l)/3 + l;

		double m2 = (r-l)/3*2 + l;

		double res1 = cal3(m1), res2 = cal3(m2);

		res = min(res, min(res1, res2));

		if(res1<res2) r = m2;

		else l = m1;

	}

	printf("%.10lf\n", res);

	return 0;

}

 

模拟退火写法，对于三维复杂度更低：

#include<cstdio>

#include<cmath>

#include<algorithm>

using namespace std;

const double eps = 1e-5;

struct Point{

    double x, y, z;

}p[110], op;

int n;

inline double dist(Point &a, Point &b) {

    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)+(a.z-b.z)*(a.z-b.z));

}

void solve() {

    double ans, delta = 10000.0;

    double maxDis, tempDis;

    while(delta>eps){

        int id = 0;

        maxDis = dist(op, p[id]);

        for(int i=1;i<n;i++){

            tempDis=dist(op,p[i]);

            if(tempDis>maxDis){

                maxDis = tempDis;

                id = i;

            }

        }

        ans = maxDis;

        op.x += (p[id].x-op.x)/maxDis*delta;

        op.y += (p[id].y-op.y)/maxDis*delta;

        op.z += (p[id].z-op.z)/maxDis*delta;

        delta *= 0.98;

    }

    printf("%.10lf\n", ans);

}

int main() {

    while(scanf("%d", &n)!=EOF && n) {

        op.x = op.y = op.z = 0;

        for(int i=0;i<n;i++) {

            scanf("%lf %lf %lf", &p[i].x, &p[i].y, &p[i].z);

        }

        solve();

    }

    return 0;

}

---

POJ2069 Super Star

这一题三分做法会T，只能用模拟退火才能过。

注意初始点选择。

#include<cstdio>

#include<cmath>

#include<algorithm>

using namespace std;

const double eps = 1e-5;

struct Point{

    double x, y, z;

}p[35], op;

int n;

inline double dist(Point &a, Point &b) {

    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)+(a.z-b.z)*(a.z-b.z));

}

void solve() {

    double ans, delta = 100.0;

    double maxDis, tempDis;

    while(delta>eps){

        int id = 0;

        maxDis = dist(op, p[id]);

        for(int i=1;i<n;i++){

            tempDis=dist(op,p[i]);

            if(tempDis>maxDis){

                maxDis = tempDis;

                id = i;

            }

        }

        ans = maxDis;

        op.x += (p[id].x-op.x)/maxDis*delta;

        op.y += (p[id].y-op.y)/maxDis*delta;

        op.z += (p[id].z-op.z)/maxDis*delta;

        delta *= 0.98;

    }

    printf("%.5lf\n", ans);

}

int main() {

    while(scanf("%d", &n)!=EOF && n) {

        op.x = op.y = op.z = 0;

        for(int i=0;i<n;i++) {

            scanf("%lf %lf %lf", &p[i].x, &p[i].y, &p[i].z);

            op.x += p[i].x;

            op.y += p[i].y;

            op.z += p[i].z;

        }

        op.x /= n; op.y /= n; op.z /= n;

        solve();

    }

    return 0;

}

---

HDU3007 HDU3932 类似。

注意HDU3932 n==1采用模拟退火要特判。。。。

D.Country Meow 最小球覆盖 三分套三分套三分 && 模拟退火的相关教程结束。