AOJ 2423 : Code Art Online

問題文

http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=2423

解法

ある人が穴を通れるかどうかは、穴の大きさが、その人を表す多角形全てを内包する円、すなわちその点集合の最小包含円以上の経をもつかどうかで判定できます。
各多角形について最小包含円を求めてから、辞書順最小になるように通る穴を決めていきます。

最小包含円を求める

全ての三点の組み合わせについて、それらによって構成される三角形の最小包含円の内、経が最大のものが全体の最小包含円となります。
三角形の最小包含円はその最大の角の大きさによって分けられ、

• 最長辺を直径とする円（鈍角三角形の場合）
• 三角形の外接円（ otherwise ）

となります。
直角三角形の場合上記二つの円は一致するので、どちらで求めてもよいでしょう。

鈍角三角形かどうかの判定ですが、鈍角の余弦は負になるので、そのような角度をもつ角が存在するかどうかで判定できます。
さらに、余弦が負であることと内積が負であることは同値です。
従って、三つの角について、その挟む辺の内積を求めることで判定できます。

外接円を求めるには、まず余弦定理 を用い、 の値を得ます。
次に を求めます。
を求めるには、 を変形して としてもよいですが、逆余弦関数 を用い、 として求めることもできます。
その後、正弦定理より として外接円の半径 r を得ます。

また、凹んでいる部分（＝凸包に含まれない点）が最小包含円に接することはないので、予め凸包を求め、凸包上の点だけを考えることで計算量を落とすことができます。

辞書順最小の割り当てを求める

辞書順最小になるよう人に穴を割り振るには、残りの「塞がってしまっても残りの人が入る穴が残る」ような穴のうちで最も番号が若いものを貪欲に割り振っていきます。
「残りの人が入れるかどうか」の判定では、最適に割り当てた場合にどうなるかだけが分かればよいので、辞書順はどうでもよく、人・穴の経だけが重要になります。
残りの人・穴をそれぞれ昇順にソートし、それぞれの人に、残りの入れる穴の内で最も小さいものを割り当てていきます。
そのような穴が存在しない場合、残りの人に穴を割り振ることはできません。

コード

typedef vector<int> VI;

#define REP( i, m, n ) for ( int i = (int)( m ); i < (int)( n ); ++i )
#define EACH( v, c ) for ( auto &v : c )

#define ALL( c ) (c).begin(), (c).end()

#define PB( n ) push_back( n )
#define EXIST( c, e ) ( (c).find( e ) != (c).end() )

const double EPS = 1e-9;

struct Point
{
	double x, y;

	Point() : x( 0 ), y( 0 ) {}
	Point( const double x, const double y ) : x( x ), y( y ) {}
	Point( const Point &a ) : x( a.x ), y( a.y ) {}
	
	Point operator + ( const Point &a ) const
	{
		return Point( x + a.x, y + a.y );
	}

	Point& operator += ( const Point &a ) 
	{
		x += a.x;
		y += a.y;

		return *this;
	}

	Point operator - ( const Point &a ) const
	{
		return Point( x + ( -a.x ), y + ( -a.y ) );
	}

	Point& operator -= ( const Point &a )
	{
		x -= a.x;
		y -= a.y;

		return *this;
	}

	Point operator * ( const double &a ) const
	{
		return Point( x * a, y * a );
	}

	Point& operator *= ( const double &a )
	{
		x *= a;
		y *= a;

		return *this;
	}

	Point operator / ( const double &a ) const
	{
		return Point( x / a, y / a );
	}

	Point& operator /= ( const double &a )
	{
		x /= a;
		y /= a;

		return *this;
	}

	bool operator < ( const Point &a ) const
	{
		return x == a.x ? y < a.y : x < a.x;
	}

	double abs() const
	{
		return sqrt( pow( x, 2. ) + pow( y, 2. ) );
	}
};

// 二点間の距離
double distance( const Point &a, const Point &b )
{
	return sqrt( pow( a.x - b.x, 2. ) + pow( a.y - b.y, 2. ) );
}

// 内積（ドット積）
double dot( const Point &a, const Point &b )
{
	return a.x * b.x + a.y * b.y;
}

// 外積（クロス積）
double cross( const Point &a, const Point &b )
{
	return a.x * b.y - a.y * b.x;
}

// 凸包を求める（ Graham Scan ）
vector< Point > convexHull( vector< Point > ps )
{
	const int N = ps.size();

	sort( ALL( ps ) );

	int k = 0;
	vector< Point > convex( N * 2 );

	for ( int i = 0; i < N; i++ )
	{
		while ( 2 <= k && cross( convex[ k - 1 ] - convex[ k - 2 ], ps[i] - convex[ k - 1 ] ) <= EPS )
		{
			k--;
		}

		convex[ k++ ] = ps[i];
	}

	for ( int i = N - 2, t = k; 0 <= i; i-- )
	{
		while ( t < k && cross( convex[ k - 1 ] - convex[ k - 2 ], ps[i] - convex[ k - 1 ] ) <= EPS )
		{
			k--;
		}

		convex[ k++ ] = ps[i];
	}

	convex.resize( k - 1 );

	return convex;
}

// 三角形の最小包含円を求める
double enclosing( Point A, Point B, Point C )
{
	double a = distance( B, C ), b = distance( C, A ), c = distance( A, B );

	if ( dot( B - A, C - A ) <= EPS ||
		 dot( A - B, C - B ) <= EPS ||
		 dot( A - C, B - C ) <= EPS ) // 内積が 0 以下 ⇔ 直角または鈍角の角がある
	{
		return max( a, max( b, c ) ) / 2;
	}
	else
	{
		double cosA = ( pow( b, 2. ) + pow( c, 2. ) - pow( a, 2. ) ) / ( 2 * b * c );
		double sinA = sin( acos( cosA ) );

		return abs( a / sinA ) / 2;
	}
}

// 各多角形について、凸包を求めてから最小包含円を求める
vector<double> convexHullRs( vector< vector< Point > > polygons )
{
	vector<double> res;

	EACH( polygon, polygons )
	{
		vector< Point > convex( convexHull( polygon ) );
		const int N = convex.size();

		double r = 0;

		REP( i, 0, N )
		{
			REP( j, i + 1, N )
			{
				REP( k, j + 1, N )
				{
					r = max( r, enclosing( convex[i], convex[j], convex[k] ) );
				}
			}
		}

		res.PB( r );
	}

	return res;
}

// 補集合を求める
template < typename T >
vector< T > diffSet( const vector< T > &whole, const vector<bool> &used )
{
	vector< T > res;

	REP( i, 0, whole.size() )
	{
		if ( used[i] )
		{
			continue;
		}

		res.PB( whole[i] );
	}

	return res;
}

// 残りの人が全員入れるか
bool valid( VI holes, vector<double> rs )
{
	sort( ALL( holes ) );
	sort( ALL( rs ) );

	int last = -1;

	REP( i, 0, rs.size() )
	{
		bool ok = false;

		REP( j, last + 1, holes.size() )
		{
			if ( rs[i] <= holes[j] + EPS )
			{
				ok = true;
				last = j;

				break;
			}
		}

		if ( !ok )
		{
			return false;
		}
	}

	return true;
}

// 辞書順最小になるように割り当てる
VI matching( const VI &holes, const vector<double> &rs )
{
	const int N = holes.size(), M = rs.size();

	VI res;
	vector<bool> closed( N, false ), entered( M, false );

	REP( i, 0, M )
	{
		REP( j, 0, N )
		{
			if ( closed[j] || entered[i] || holes[j] < rs[i] - EPS )
			{
				continue;
			}

			vector<bool> nextClosed( closed ), nextEntered( entered );
			nextClosed[j] = true;
			nextEntered[i] = true;

			if ( valid( diffSet( holes, nextClosed ), diffSet( rs, nextEntered ) ) )
			{
				entered = nextEntered;
				closed = nextClosed;
				res.PB( j + 1 );

				break;
			}
		}
	}

	return res.size() == M ? res : VI();
}

int main()
{
	cin.tie( 0 );
	ios::sync_with_stdio( false );

	int n, m;
	cin >> n >> m;

	VI holes( n );

	EACH( hole, holes )
	{
		cin >> hole;
	}

	vector< vector< Point > > players( m );

	EACH( player, players )
	{
		int p;
		cin >> p;

		player.resize( p );

		EACH( point, player )
		{
			cin >> point.x >> point.y;
		}
	}

	VI res( matching( holes, convexHullRs( players ) ) );

	if ( res.size() < m )
	{
		cout << "NG" << endl;
	}
	else
	{
		EACH( tmp, res )
		{
			cout << tmp << endl;
		}
	}

	return 0;
}