HDU 3644 计算几何 模拟退火

/*
 * =====================================================================================
 *
 *       Filename:  hdu3644.cpp
 *
 *    Description:  Calculate Geometry -- Simulated annealing
 *
 *        Version:  1.0
 *        Created:  2012年01月21日 10时41分15秒
 *       Revision:  none
 *       Compiler:  gcc
 *
 *         Author:  SphinX (Yuxiang Ye), sphinx2010@yahoo.cn
 *   Organization:  
 *
 * =====================================================================================
 */

/*
*    10年杭州网络赛的一道题，碉堡了～貌似是我A的最艰难的题了，虽然一开始
*就看出是模拟退火，但无奈计算几何基础太差，敲的过程中各种错误，然后又重新学
*了遍计算几何基本功，但还是各种WA跟TLE，最后发现是精度问题。。。
*    其中辛苦一言难尽。。。
*/

#include	<stdlib.h>
#include	<iostream>
#include	<algorithm>
#include	<cstdio>
#include	<cmath>
using namespace std;

typedef double LF;

const LF eps = 1e-3;

const int MAXN = 104;

struct Point 
{
	LF x, y;
	Point(){}
	Point(const LF _x, const LF _y)
	{
		x = _x;
		y = _y;
	}
};				/* ----------  end of struct Point  ---------- */

typedef struct Point Point;

struct Segment 
{
	Point ps, pe;
	Segment(){}
	Segment(const Point _ps, const Point _pe)
	{
		ps = _ps;
		pe = _pe;
	}
};				/* ----------  end of struct Segment  ---------- */

typedef struct Segment Segment;

int n;
LF radius, min_dist[MAXN];
Point pnt[MAXN], tmp[MAXN];
Segment seg[MAXN];


/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  dblcmp
 *  Description:  
 * =====================================================================================
 */
	int
dblcmp ( const LF & x )
{
	return fabs(x) < eps ? 0 : (x < 0.0 ? -1 : 1);
}		/* -----  end of function dblcmp  ----- */

/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  xmult
 *  Description:  
 * =====================================================================================
 */
	LF
xmult ( const Point & p1, const Point & p2, const Point & p0 )
{
	LF vx[2] = {p1.x - p0.x, p2.x - p0.x};
	LF vy[2] = {p1.y - p0.y, p2.y - p0.y};
	return vx[0] * vy[1] - vx[1] * vy[0];
}		/* -----  end of function xmult  ----- */

/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  pnt2pnt_dist
 *  Description:  
 * =====================================================================================
 */
	LF
pnt2pnt_dist ( const Point & pa, const Point &pb )
{
	return sqrt((pa.x - pb.x) * (pa.x - pb.x)
		+ (pa.y - pb.y) * (pa.y - pb.y));
}		/* -----  end of function pnt2pnt_dist  ----- */

/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  pnt2seg_dist
 *  Description:  
 * =====================================================================================
 */
	LF
pnt2seg_dist ( const Point & pt, const Segment & seg )
{
	LF a = pnt2pnt_dist(pt, seg.pe);
	if( a < eps )
	{
		return a;
	}
	LF b = pnt2pnt_dist(pt, seg.ps);
	if( b < eps )
	{
		return b;
	}
	LF c = pnt2pnt_dist(seg.pe, seg.ps);
	if( c < eps )
	{
		return a;
	}
	if( a * a >= b * b + c * c )
	{
		return b;
	}
	if( b * b >= a * a + c * c )
	{
		return a;
	}
	LF l = (a + b + c) * 0.5;
	LF s = sqrt(l * (l - a) * (l - b) * (l - c));
	return 2.0 * s / c;
}		/* -----  end of function pnt2seg_dist  ----- */

/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  pnt2poly_dist
 *  Description:  
 * =====================================================================================
 */
	LF
pnt2poly_dist ( const Point & pt )
{
	LF mn = pnt2seg_dist(pt, seg[0]);
	for( int i = 1; i < n ; ++i )
	{
		mn = min(mn, pnt2seg_dist(pt, seg[i]));
	}
	
	return mn;
}		/* -----  end of function pnt2poly_dist  ----- */

/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  seg_cross_judge
 *  Description:  
 * =====================================================================================
 */
	bool
seg_cross_judge ( const Segment & u, const Segment & v )
{
	return max(u.pe.x, u.ps.x) >= min(v.pe.x, v.ps.x) 
			&& max(v.pe.x, v.ps.x) >= min(u.pe.x, u.ps.x)
			&& max(u.pe.y, u.ps.y) >= min(v.pe.y, v.ps.y)
			&& max(v.pe.y, v.ps.y) >= min(u.pe.y, u.ps.y)
			&& -1 == dblcmp(xmult(v.ps, u.pe, u.ps)) * dblcmp(xmult(v.pe, u.pe, u.ps))
			&& -1 == dblcmp(xmult(u.ps, v.pe, v.ps)) * dblcmp(xmult(u.pe, v.pe, v.ps));
}		/* -----  end of function seg_cross_judge  ----- */


/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  pnt_inside_judge
 *  Description:  
 * =====================================================================================
 */
	bool
pnt_inside_judge ( const Point & pt )
{
	LF ang = rand() + 0.0;
	Point ps = Point(pt.x + 10000.0 * cos(ang), pt.y + 10000.0 * sin(ang));
	Segment sg = Segment(pt, ps);
	int cnt = 0;
	for( int i = 0; i < n ; ++i )
	{
		if (seg_cross_judge(sg, seg[i]))
		{
			++cnt;
		}
	}
	return cnt % 2;
}		/* -----  end of function pnt_inside_judge  ----- */

/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  exist_judge
 *  Description:  模拟退火
 * =====================================================================================
 */
	bool
exsist_judge ()
{
	LF maxr = 0.0;
	for( int i = 0; i < n ; ++i )
	{
		tmp[i].x = (pnt[i].x + pnt[i + 1].x) * 0.5;
		tmp[i].y = (pnt[i].y + pnt[i + 1].y) * 0.5;
		min_dist[i] = pnt2poly_dist(tmp[i]);
	}
	for (int i = 0; i < n - 1; ++i)
	{
		for (int j = i + 1; j < n; ++j)
		{
			maxr = max(maxr, pnt2pnt_dist(pnt[i], pnt[j]));
		}
	}
	while (maxr > 1e-3)
	{
		for( int k = 0; k < 5 ; ++k )
		{
			for (int i = 0; i < n; ++i)
			{
				LF ang = rand() + 0.0;
				Point pk = Point(tmp[i].x + maxr * cos(ang), tmp[i].y + maxr * sin(ang));
				if (!pnt_inside_judge(pk))
				{
					continue ;
				}
				LF buf = pnt2poly_dist(pk);
				if (buf + eps > radius)
				{
					return true;
				}
				if (buf > min_dist[i])
				{
					min_dist[i] = buf;
					tmp[i] = pk;
				}
			}
		}
		maxr *= 0.9;
	}
	
	return false;
}		/* -----  end of function exist_judge  ----- */
/* 
 * ===  FUNCTION  ======================================================================
 *         Name:  main
 *  Description:  VIM自动生成的代码，比较吓人～
 * =====================================================================================
 */
	int
main ( int argc, char *argv[] )
{
	srand(static_cast<unsigned int>(time(NULL)));
	while( scanf ( "%d", &n ), n )
		{
			for( int i = 0; i < n ; ++i )
			{
				scanf ( "%lf %lf", &pnt[i].x, &pnt[i].y );
			}
			scanf ( "%lf", &radius );
			pnt[n] = pnt[0];
			for( int i = 0; i < n ; ++i )
			{
				seg[i] = Segment(pnt[i], pnt[i + 1]);
			}
			cout << (exsist_judge() ? "Yes" : "No") << endl;
		}
	return EXIT_SUCCESS;
}				/* ----------  end of function main  ---------- */