問題E『六角形』 : MUJIN プログラミングチャレンジ

ほとんど数学の問題でした。

問題文 : 六角形

概要

上の凸六角形(頂点は反時計回り)が与えられる。
この六角形を断面として持つ直方体は存在するか。
存在しなければ-1を答え、存在すればそのような直方体の体積の最小値を答えよ。

制約

座標はいずれも整数。

解説

実は、存在するときは直方体は一意に定まる。

直方体の断面が六角形になるのは、図のような場合。

まず、六角形の対辺は平行でなければならない。
すなわち、P[6]=P[0]として、P[i]P[i+1]//P[i+3]P[i+4]がi=0,1,2について成立しなければならない。
そうでなければ-1。

次に面を決めていく。まず、P[2]P[3]を通る面を考える。
面は、通る点と法線ベクトルによって決まるので、あとは単位法線ベクトルn[23]を与えればよい。
線分P[2]P[3]はこの面の上にあるので、ベクトルP[2]P[3]はn[23]と直交する。
つまり、でなければならない。

これを満たすn[23]は、P[2]P[3]を軸として、回転の分だけ自由度がある。
xy平面とn[23]のなす角をθとする。xy平面に対する鏡映対称性を用いて、としてよい。
つまり、P[2]P[3]と直交するxy平面上の単位ベクトルをv23nと名付け、z軸方向の単位ベクトルをezとすると、と書ける。

他の面はどうか。隣の、P[4]P[5]を通る面を考えよう。
やはり単位法線ベクトルn[45]を決めればよいが、今度は一意に定まってしまう。
なぜなら、n[45]はP[4]P[5]だけでなく、n[23]とも直交しなければならないからだ。
ゆえに、である。

同様に、P[0]P[1]を通る面も決まる。単位法線ベクトルn[01]は、である。

これで直方体のすべての面が決まった。最後に、これがちゃんと直方体をなしていることを確認しなければならない。
つまり、n[45]とn[01]が直交しているかを調べる必要がある。
そのために、を解く必要がある。
公式を用いて計算すると、
となる。

この右辺をKとして、でなければ直方体は作れないので、-1となる。
そうでなければ、θが定まり、すべての法線ベクトルが定まる。
よって、あとは各頂点を求めれば、体積が求まる。
最小値も何も、存在するときは一意に定まってしまうのである。

例えば、図の点Hは、P[3]から、P[4]P[5]を通る面に下した垂線の足である。
よって、H=P[3]+x*n[45]と置けて、P[4]Hがn[45]と直交することから、xが求まり、
となる。

同様に、点G,C,Bも求まるので、あとは各辺の長さをかければよい。

ソースコード

using System;
using System.Linq;
class K { static void Main() { new E(); } }
class E
{
	int[] G() { return Console.ReadLine().Split().Select(_ => int.Parse(_)).ToArray(); }
	Vector[] P;
	public E()
	{
		P = new Vector[7];
		for (var i = 0; i < 6; i++) { var I = G(); P[i] = new Vector(I[0], I[1]); }
		P[6] = P[0];
		Console.WriteLine(Enumerable.Range(0, 3).All(i => IsParallel(Vec(i, i + 1), Vec(i + 3, i + 4))) ? Volume() : -1);
	}
	Vector Vec(int from, int to) { return P[to] - P[from]; }
	bool IsParallel(Vector a, Vector b) { return (a % b).Abs == 0; }
	Vector Project(Vector a, Vector b, Vector n) { return a + ((b - a) * n) * n; }
	double Volume()
	{
		var ez = new Vector(0, 0, 1);
		// % は外積
		var v23n = (Vec(2, 3) % ez).Unit;
		var v45 = Vec(4, 5);
		var v01 = Vec(0, 1);

		// n01 と n45 が直交する角度を求める
		var K = (v45 * v01) / ((v45 * v23n) * (v01 * v23n));
		if (K <= 0 || 1 <= K) return -1;
		var c = Math.Sqrt(K);
		var s = Math.Sqrt(1 - K);

		// 面の法線ベクトル
		var n23 = c * v23n + s * ez;
		var n45 = (v45 % n23).Unit;
		var n01 = (v01 % n23).Unit;

		// 直方体の頂点
		var H = Project(P[3], P[4], n45);
		var G = Project(P[3], P[2], n45);
		var C = Project(P[2], P[1], n01);
		var B = Project(P[1], P[0], n23);
		return (H - G).Abs * (G - C).Abs * (C - B).Abs;
	}
}
class Vector
{
	public double X, Y, Z;
	public Vector(double x = 0, double y = 0, double z = 0) { X = x; Y = y; Z = z; }
	public double Abs { get { return Math.Sqrt(this * this); } }
	public Vector Unit { get { return (1 / Abs) * this; } }
	public static Vector operator -(Vector p, Vector q) { return new Vector(p.X - q.X, p.Y - q.Y, p.Z - q.Z); }
	public static Vector operator +(Vector p, Vector q) { return new Vector(p.X + q.X, p.Y + q.Y, p.Z + q.Z); }
	public static Vector operator *(double r, Vector p) { return new Vector(p.X * r, p.Y * r, p.Z * r); }
	public static double operator *(Vector p, Vector q) { return p.X * q.X + p.Y * q.Y + p.Z * q.Z; }
	// 外積
	public static Vector operator %(Vector p, Vector q) { return new Vector(p.Y * q.Z - p.Z * q.Y, p.Z * q.X - p.X * q.Z, p.X * q.Y - p.Y * q.X); }
}

感想

ただの数学でした。
最初はn[01]とn[45]が勝手に直交してくれると思い込んで詰まり、
それに気づいてからも、なぜか直交性を判定するのにで計算し、体積が0になってしまい…。
結局、終了8分後にAC。