fork download

copy

namespace h3d {
    using FLOAT = double;
	template <unsigned int N> class Matrix;
	template <unsigned int N> class Vector;
};
 
namespace h3d {
	/**
	 * 拡大・縮小変換行列を作る。
	 * @param coef 各座標の係数(COEFficient)。
	 * coef > 1.0 なら拡大, coef < 1.0 なら縮小。
	 * @return 作った拡大・縮小変換行列。
	 */
	Matrix<3> scaling(const Vector<3>& coef);
 
	/**
	 * 回転変換行列を作る。
	 * @param axis この軸の周りを回転する。
	 * @param rad 回転する弧度(RADian)。
	 * @return 作った回転変換行列。
	 */
	Matrix<3> rotation(const Vector<3>& axis, const FLOAT& rad);
 
	/**
	 * 平行移動変換行列を作る。
	 * @param os 変換後の位置から変換前の位置を引いた差分(OffSet)。
	 * @return 作った平行移動変換行列。
	 */
	Matrix<4> translation(const Vector<3>& os);
};
 
#include <iostream>
 
namespace h3d {
	using std::ostream;
	using std::string;
 
	/**
	 * N次元正方行列クラス。
	 * @param N 行列の次元数。
	 */
	template <unsigned int N> 
	class Matrix {
	public:
		/** 成分(Component)の配列。N×N個。 */
		FLOAT c[N][N];
 
		/** 単位行列。 */
		static const Matrix IDENTITY;
 
		/**
		 * 各成分をゼロクリアする。
		 */
		inline void clear();
 
		/**
		 * 行列式(DETerminant)を計算する。
		 * @return 計算した行列式。
		 */
		inline FLOAT det() const;
 
		/**
		* 次元を一つ拡張(EXTend)する。
		 * @param mat この行列から最後の行と列の成分をコピーする。
		 * @return 拡張した行列。
		 */
		inline Matrix<N + 1> ext(const Matrix<N + 1>& mat = Matrix<N + 1>::IDENTITY) const;
 
		/**
		 * 逆行列(INVerse matrix)を計算する。
		 * @return 計算した逆行列。
		 * @throw string 行列が正則ではないならメッセージをスローする。
		 */
		inline Matrix inv() const throw(string);
 
		/**
		 * 二つの行列が等しくないかどうかを判定する。
		 * @param rhs 右側の行列。
		 * @return 等しくないならtrue, 等しいならfalse。
		 */
		inline bool operator !=(const Matrix& rhs) const;
 
		/**
		 * 行列で乗算する。
		 * @param rhs 乗数となる右側の行列。
		 * @return 乗算した行列。
		 */
		inline Matrix operator *(const Matrix& rhs) const;
 
		/**
		 * ベクトルとの積を計算する。
		 * @param rhs 乗数となる右側のベクトル。
		 * @return 計算したベクトル。
		 */
		inline Vector<N> operator *(const Vector<N>& rhs) const;
 
		/**
		 * スカラー値で乗算する。
		 * @param rhs 乗数となる右側のスカラー値。
		 * @return 乗算した行列。
		 */
		inline Matrix operator *(const FLOAT& rhs) const;
 
		/**
		 * 行列で乗算して、代入する。
		 * @param rhs 乗数となる右側の行列。
		 * @return 乗算して、代入した行列。
		 */
		inline Matrix operator *=(const Matrix& rhs);
 
		/**
		 * スカラー値で乗算して、代入する。
		 * @param rhs 乗数となる右側のスカラー値。
		 * @return 乗算して、代入した行列。
		 */
		inline Matrix operator *=(const FLOAT& rhs);
 
		/**
		 * 行列を加算する。
		 * @param rhs 加数となる右側の行列。
		 * @return 加算した行列。
		 */
		inline Matrix operator +(const Matrix& rhs) const;
 
		/**
		 * 行列を加算して、代入する。
		 * @param rhs 加数となる右側の行列。
		 * @return 加算して、代入した行列。
		 */
		inline Matrix operator +=(const Matrix& rhs);
 
		/**
		 * すべての成分の符号を反転する。
		 * @return 符号を反転した行列。
		 */
		inline Matrix operator -() const;
 
		/**
		 * 行列を減算する。
		 * @param rhs 減数となる右側の行列。
		 * @return 減算した行列。
		 */
		inline Matrix operator -(const Matrix& rhs) const;
 
		/**
		 * 行列を減算して、代入する。
		 * @param rhs 減数となる右側の行列。
		 * @return 減算して、代入した行列。
		 */
		inline Matrix operator -=(const Matrix& rhs);
 
		/**
		 * スカラー値で除算する。
		 * @param rhs 除数となる右側のスカラー値。
		 * @return 除算した行列。
		 */
		inline Matrix operator /(const FLOAT& rhs) const;
 
		/**
		 * スカラー値で除算して、代入する。
		 * @param rhs 除数となる右側のスカラー値。
		 * @return 除算して、代入した行列。
		 */
		inline Matrix operator /=(const FLOAT& rhs);
 
		/**
		 * 二つの行列が等しいかどうかを判定する。
		 * @param rhs 右側の行列。
		 * @return 等しいならtrue, 等しくないならfalse。
		 */
		inline bool operator ==(const Matrix& rhs) const;
 
		/**
		 * 行と列を取り除いて、小行列(SUBmatrix)を作る。
		 * @param row 取り除く行。
		 * @param col 取り除く列。
		 * @return 作った行列。
		 */
		inline Matrix<N - 1> sub(const unsigned int& row, const unsigned int& col) const;
 
		/**
		 * 次元を一つ切り縮める(TRUNCate)。
		 * @return 切り縮めた行列。
		 */
		inline Matrix<N - 1> trunc() const;
	};
 
	/**
	 * 行列を出力する。
	 * @param N 行列の次元数。
	 * @param out 出力先となるストリーム。
	 * @param mat 出力する行列。
	 * @return 出力したストリーム。
	 */
	template <unsigned int N> 
	inline ostream& operator <<(ostream& out, const Matrix<N>& mat);
};
 
#include <iostream>
 
namespace h3d {
	using std::copy;
	using std::ostream;
 
	/**
	 * N次元ベクトルクラス。
	 * @param N ベクトルの次元数。
	 */
	template <unsigned int N> 
	class Vector {
	public:
		/** 成分(Component)の配列。N個。 */
		FLOAT c[N];
 
		/**
		 * 各成分をゼロクリアする。
		 */
		inline void clear();
 
		/**
		 * 次元を一つ拡張(EXTend)する。
		 * @param c 追加する成分。
		 * @return 拡張したベクトル。
		 */
		inline Vector<N + 1> ext(const FLOAT& c = 1.0) const;
 
		/**
		 * 内積(INner PROduct)を計算する。
		 * @param rhs 右側のベクトル。
		 * @return 計算した内積。
		 */
		inline FLOAT inpro(const Vector& rhs) const;
 
		/**
		 * 二つのベクトルが等しくないかどうかを判定する。
		 * @param rhs 右側のベクトル。
		 * @return 等しくないならtrue, 等しいならfalse。
		 */
		inline bool operator !=(const Vector& rhs) const;
 
		/**
		 * スカラー値で乗算する。
		* @param rhs 乗数となる右側のスカラー値。
		 * @return 乗算したベクトル。
		 */
		inline Vector operator *(const FLOAT& rhs) const;
 
		/**
		 * スカラー値で乗算して、代入する。
		 * @param rhs 乗数となる右側のスカラー値。
		 * @return 乗算して、代入したベクトル。
		 */
		inline Vector& operator *=(const FLOAT& rhs);
 
		/**
		 * ベクトルを加算する。
		 * @param rhs 加数となる右側のベクトル。
		 * @return 加算したベクトル。
		 */
		inline Vector operator +(const Vector& rhs) const;
 
		/**
		 * ベクトルを加算して、代入する。
		 * @param rhs 加数となる右側のベクトル。
		 * @return 加算して、代入したベクトル。
		 */
		inline Vector& operator +=(const Vector& rhs);
 
		/**
		 * 各成分の符号を反転する。
		 * @return 符号を反転したベクトル。
		 */
		inline Vector operator -() const;
 
		/**
		 * ベクトルを減算する。
		 * @param rhs 減数となる右側のベクトル。
		 * @return 減算したベクトル。
		 */
		inline Vector operator -(const Vector& rhs) const;
 
		/**
		 * ベクトルを減算して、代入する。
		 * @param rhs 減数となる右側のベクトル。
		 * @return 減算して、代入したベクトル。
		 */
		inline Vector& operator -=(const Vector& rhs);
 
		/**
		 * スカラー値で除算する。
		 * @param rhs 除数となる右側のスカラー値。
		 * @return 除算したベクトル。
		 */
		inline Vector operator /(const FLOAT& rhs) const;
 
		/**
		 * スカラー値で除算して、代入する。
		 * @param rhs 除数となる右側のスカラー値。
		 * @return 除算して、代入したベクトル。
		 */
		inline Vector& operator /=(const FLOAT& rhs);
 
		/**
		 * 二つのベクトルが等しいかどうかを判定する。
		 * @param rhs 右側のベクトル
		 * @return 等しいならtrue, 等しくないならfalse。
		 */
		inline bool operator ==(const Vector& rhs) const;
 
		/**
		 * 3次元ベクトルの外積(OUTer PROduct)を計算する。
		 * @param rhs 右側のベクトル。
		 * @return 計算した外積。
		 */
		inline Vector outpro(const Vector& rhs) const;
 
		/**
		 * 次元を一つ切り縮める(TRUNCate)。
		 * @return 切り縮めたベクトル。
		 */
		inline Vector<N - 1> trunc() const;
	};
 
	/**
	 * ベクトルを出力する。
	 * @param N ベクトルの次元数。
	 * @param out 出力先となるストリーム。
	 * @param vec 出力するベクトル。
	 * @return 出力したストリーム。
	 */
	template <unsigned int N> 
	inline ostream& operator <<(ostream& out, const Vector<N>& vec);
};
 
namespace h3d {
	/**
	 * 角度(ANGle)を弧度(RADian)に変換する。
	 * @param ang 変換する角度(ANGle)。
	 * @return 変換した弧度。
	 */
	inline FLOAT ang2rad(const FLOAT& ang);
 
	/**
	 * 二つの数が近いかどうかを判定する。
	 * @param x 一つ目の数。
	 * @param y 二つ目の数。
	 * @return 近いならtrue, 遠いならfalse。
	 */
	inline bool nearEqual(const FLOAT& x, const FLOAT& y);
};
 
#include <iostream>
#include <string>
 
namespace h3d {
	using std::cout;
	using std::endl;
	using std::string;
 
	#define ASSERT(pred) assert(#pred, (pred));
 
	#define PRINT(val) cout << #val << "=" << (val) << endl;
 
	/**
	 * アサーションを実行する。
	 * 成功なら標準出力に結果を出力する。
	 * @param pred_str 判定する式を表した文字列。
	 * @param pred_res 判定結果。trueなら成功，falseなら失敗と判定する。
	 * @throw string 失敗ならメッセージをスローする。
	 */
	inline void assert(const string& pred_str, const bool& pred_res) throw(string);
};
 
namespace h3d {
	class Test { public: virtual void run() const = 0; };
};
 
#include <memory>
#include <vector>
 
namespace h3d {
	using std::shared_ptr;
	using std::vector;
 
	class TestSet {
	public:
		TestSet();
		void run() const;
	protected:
		vector<shared_ptr<Test>> tests;
	};
};
 
namespace h3d {
	class Test1 : public Test {
	public:
		virtual void run() const override;
	};
};
 
namespace h3d {
	class Test2 : public Test {
	public:
		virtual void run() const override;
	};
};
 
namespace h3d {
	class Test3 : public Test {
	public:
		virtual void run() const override;
	};
};
 
////////////////////////////////////////////////////////////////////////////////
 
#include <cmath>
 
namespace h3d {
	Matrix<3> scaling(const Vector<3>& coef) {
		Matrix<3> res(Matrix<3>::IDENTITY);
		for (auto i = 0; i < 3; i++) res.c[i][i] = coef.c[i];
		return res;
	}
 
	Matrix<3> rotation(const Vector<3>& axis, const FLOAT& rad) {
		FLOAT xx = axis.c[0] * axis.c[0];
		FLOAT xy = axis.c[0] * axis.c[1];
		FLOAT xz = axis.c[0] * axis.c[2];
		FLOAT yy = axis.c[1] * axis.c[1];
		FLOAT yz = axis.c[1] * axis.c[2];
		FLOAT zz = axis.c[2] * axis.c[2];
		FLOAT s = sin(rad);
		FLOAT c = cos(rad);
		FLOAT d = 1.0 - c;
		FLOAT xxd = xx * d;
		FLOAT xyd = xy * d;
		FLOAT zs = axis.c[2] * s;
		FLOAT xzd = xz * d;
		FLOAT ys = axis.c[1] * s;
		FLOAT yyd = yy * d;
		FLOAT yzd = yz * d;
		FLOAT xs = axis.c[0] * s;
		FLOAT zzd = zz * d;
		return Matrix<3>{
			c + xxd,	xyd - zs,	xzd + ys, 
			xyd + zs,	c + yyd,	yzd - xs, 
			xzd - ys,	yzd + xs,	c + zzd, 
		};
	}
 
	Matrix<4> translation(const Vector<3>& os) {
		Matrix<4> res(Matrix<4>::IDENTITY);
		for (auto i = 0; i < 3; i++) res.c[i][3] = os.c[i];
		return res;
	}
};
 
namespace h3d {
	template <unsigned int N> 
	inline void Matrix<N>::clear() {
		// 各成分をゼロクリアする。
		for (auto row = 0; row < N; row++) {
			for (auto col = 0; col < N; col++) 
				this->c[row][col] = 0.0;
		}
	}
 
	template <unsigned int N> 
	inline FLOAT Matrix<N>::det() const {
		FLOAT res = 0.0;
		// 1行目の各成分ごとに反復する。
		for (auto col = 0; col < N; col++) {
			// 成分に、それと対応する小行列式を乗算する。
			FLOAT cofac = c[0][col] * sub(0, col).det();
			// (-1)^col
			if ((col & 0x1) == 1) cofac = -cofac;
			// 結果に余因子を加算する。
			res += cofac;
		}
		return res;
	}
 
	template <> 
	inline FLOAT Matrix<2>::det() const {
		// 再帰的な行列式計算の終着点。
		// 2×2の行列式を計算する。
		return this->c[0][0] * this->c[1][1] - this->c[0][1] * this->c[1][0];
	}
 
	template <unsigned int N> 
	inline Matrix<N + 1> Matrix<N>::ext(const Matrix<N + 1>& mat) const {
		Matrix<N + 1> res;
		// 結果の各成分ごとに反復する。
		for (auto row = 0; row < N + 1; row++) {
			for (auto col = 0; col < N + 1; col++) {
				// 行と列がN以内ならこの行列, それ以外なら引数の行列から成分をコピーする。
				if (row < N && col < N) res.c[row][col] = this->c[row][col];
				else res.c[row][col] = mat.c[row][col];
			}
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::inv() const throw(string) {
		// まず行列式を計算し、正則であることを確かめる。
		FLOAT det_res = det();
		if (nearEqual(det_res, 0.0)) throw "行列が正則ではない。";
		// 行列式の逆数を計算しておく。
		FLOAT inv_det = 1.0 / det_res;
		Matrix<N> res;
		// 各成分ごとに反復する。
		for (auto row = 0; row < N; row++) {
			for (auto col = 0; col < N; col++) {
				// 行列式の逆数に、対角の小行列式を乗算する。
				res.c[row][col] = inv_det * sub(col, row).det();
				// (-1)^(row + col)
				if (((row + col) & 0x1) == 1) 
					res.c[row][col] = -res.c[row][col];
			}
		}
		return res;
	}
 
	template <unsigned int N> 
	inline bool Matrix<N>::operator !=(const Matrix& rhs) const {
		return !operator ==(rhs);
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator *(const Matrix& rhs) const {
		Matrix<N> res;
		// 結果の各成分ごとに反復する。
		for (auto row = 0; row < N; row++) {
			for (auto col = 0; col < N; col++) {
				res.c[row][col] = 0.0;
				// 左側は対応する行の各成分, 右側は対応する列の各成分ごとに反復して、乗算する。
				for (auto i = 0; i < N; i++) 
					res.c[row][col] += this->c[row][i] * rhs.c[i][col];
			}
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Vector<N> Matrix<N>::operator *(const Vector<N>& rhs) const {
		Vector<N> res;
		// ベクトルの各成分ごとに反復する。
		for (auto i = 0; i < N; i++) {
			res.c[i] = 0.0;
			// 左側の行列は、対応する行の各成分ごとに反復して、乗算する。
			for (auto j = 0; j < N; j++) 
				res.c[i] += this->c[i][j] * rhs.c[j];
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator *(const FLOAT& rhs) const {
		Matrix<N> res;
		// 各成分にスカラー値を乗算する。
		for (auto row = 0; row < N; row++) {
			for (auto col = 0; col < N; col++) 
				res.c[row][col] = this->c[row][col] * rhs;
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator *=(const Matrix& rhs) {
		return *this = operator *(rhs);
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator *=(const FLOAT& rhs) {
		return *this = operator *(rhs);
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator +(const Matrix& rhs) const {
		Matrix<N> res;
		// 対応する成分同士で加算する。
		for (auto row = 0; row < N; row++) {
			for (auto col = 0; col < N; col++) 
				res.c[row][col] = this->c[row][col] + rhs.c[row][col];
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator +=(const Matrix& rhs) {
		return *this = operator +(rhs);
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator -() const {
		Matrix<N> res;
		// 各成分の符号を反転する。
		for (auto row = 0; row < N; row++) {
			for (auto col = 0; col < N; col++) 
				res.c[row][col] = -this->c[row][col];
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator -(const Matrix& rhs) const {
		return operator +(-rhs);
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator -=(const Matrix& rhs) {
		return *this = operator -(rhs);
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator /(const FLOAT& rhs) const {
		return operator *(1.0 / rhs);
	}
 
	template <unsigned int N> 
	inline Matrix<N> Matrix<N>::operator /=(const FLOAT& rhs) {
		return *this = operator /(rhs);
	}
 
	template <unsigned int N> 
	bool Matrix<N>::operator ==(const Matrix& rhs) const {
		bool res = true;
		if (&rhs != this) {
			// 各成分ごとに反復する。
			for (auto row = 0; row < N; row++) {
				for (auto col = 0; col < N; col++) {
					// 対応する成分同士が近ければ等しいと判定する。
					if (!nearEqual(this->c[row][col], rhs.c[row][col])) {
						res = false;
						break;
					}
				}
				if (!res) break;
			}
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Matrix<N - 1> Matrix<N>::sub(const unsigned int& row, const unsigned int& col) const {
		Matrix<N - 1> sub;
		auto sub_row = 0;
		// この行列の各成分ごとに反復する。
		// 取り除く行と列については処理をスキップする。
		for (auto sup_row = 0; sup_row < N; sup_row++) {
			if (sup_row == row) continue;
			auto sub_col = 0;
			for (auto sup_col = 0; sup_col < N; sup_col++) {
				if (sup_col == col) continue;
				// 対応する成分をコピーする。
				sub.c[sub_row][sub_col] = this->c[sup_row][sup_col];
				sub_col++;
			}
			sub_row++;
		}
		return sub;
	}
 
	template <unsigned int N> 
	inline Matrix<N - 1> Matrix<N>::trunc() const {
		Matrix<N - 1> res;
		// 結果の各成分ごとに反復する。
		for (auto row = 0; row < N - 1; row++) {
			// 対応する成分をコピーする。
			for (auto col = 0; col < N - 1; col++) 
				res.c[row][col] = this->c[row][col];
		}
		return res;
	}
 
	template <unsigned int N> 
	inline ostream& operator <<(ostream& out, const Matrix<N>& mat) {
		out << "{";
		for (auto row = 0; row < N; row++) {
			out << "{";
			for (auto col = 0; col < N; col++) out << mat.c[row][col] << ", ";
			out << "}, ";
		}
		out << "}";
	}
};
 
namespace h3d {
	template <> 
	const Matrix<2> Matrix<2>::IDENTITY{
		1.0, 0.0, 
		0.0, 1.0, 
	};
 
	template <> 
	const Matrix<3> Matrix<3>::IDENTITY{
		1.0, 0.0, 0.0, 
		0.0, 1.0, 0.0, 
		0.0, 0.0, 1.0, 
	};
 
	template <> 
	const Matrix<4> Matrix<4>::IDENTITY{
		1.0, 0.0, 0.0, 0.0, 
		0.0, 1.0, 0.0, 0.0, 
		0.0, 0.0, 1.0, 0.0, 
		0.0, 0.0, 0.0, 1.0, 
	};
};
 
#include <algorithm>
#include <iostream>
 
namespace h3d {
	using std::ostream;
 
	template <unsigned int N> 
	inline void Vector<N>::clear() {
		// 各成分をゼロクリアする。
		for (auto i = 0; i < N; i++) this->c[i] = 0.0;
	}
 
	template <unsigned int N> 
	inline Vector<N + 1> Vector<N>::ext(const FLOAT& c) const {
		Vector<N + 1> res;
		// このベクトルからN個の成分をコピーして、最後は引数の成分をコピーする。
		copy(this->c, this->c + N, res.c);
		res.c[N] = c;
		return res;
	}
 
	template <unsigned int N> 
	inline FLOAT Vector<N>::inpro(const Vector& rhs) const {
		FLOAT res = 0;
		// 対応する成分同士の積を求めて、合計する。
		for (auto i = 0; i < N; i++) res += this->c[i] * rhs.c[i];
		return res;
	}
 
	template <unsigned int N> 
	inline bool Vector<N>::operator !=(const Vector& rhs) const {
		return !operator ==(rhs);
	}
 
	template <unsigned int N> 
	inline Vector<N> Vector<N>::operator *(const FLOAT& rhs) const {
		Vector<N> res;
		// 各成分にスカラー値を乗算する。
		for (auto i = 0; i < N; i++) res.c[i] = this->c[i] * rhs;
		return res;
	}
 
	template <unsigned int N> 
	inline Vector<N>& Vector<N>::operator *=(const FLOAT& rhs) {
		return *this = operator *(rhs);
	}
 
	template <unsigned int N> 
	inline Vector<N> Vector<N>::operator +(const Vector& rhs) const {
		Vector<N> res;
		// 対応する成分同士で加算する。
		for (auto i = 0; i < N; i++) res.c[i] = this->c[i] + rhs.c[i];
		return res;
	}
 
	template <unsigned int N> 
	inline Vector<N>& Vector<N>::operator +=(const Vector& rhs) {
		return *this = operator +(rhs);
	}
 
	template <unsigned int N> 
	inline Vector<N> Vector<N>::operator -() const {
		Vector<N> res;
		// 各成分の符号を反転する。
		for (auto i = 0; i < N; i++) res.c[i] = -this->c[i];
		return res;
	}
 
	template <unsigned int N> 
	inline Vector<N> Vector<N>::operator -(const Vector& rhs) const {
		return operator +(-rhs);
	}
 
	template <unsigned int N> 
	inline Vector<N>& Vector<N>::operator -=(const Vector& rhs) {
		return *this = operator -(rhs);
	}
 
	template <unsigned int N> 
	inline Vector<N> Vector<N>::operator /(const FLOAT& rhs) const {
		return operator *(1.0 / rhs);
	}
 
	template <unsigned int N> 
	inline Vector<N>& Vector<N>::operator /=(const FLOAT& rhs) {
		return *this = operator /(rhs);
	}
 
	template <unsigned int N> 
	inline bool Vector<N>::operator ==(const Vector& rhs) const {
		bool res = true;
		if (&rhs != this) {
			// 各成分ごとに反復する。
			for (auto i = 0; i < N; i++) {
				// 対応する成分同士が近ければ等しいと判定する。
				if (!nearEqual(this->c[i], rhs.c[i])) {
					res = false;
					break;
				}
			}
		}
		return res;
	}
 
	template <unsigned int N> 
	inline Vector<N> Vector<N>::outpro(const Vector& rhs) const {
		// 各成分について、対応する列を除いた2×2の行列式を計算する。
		// P×Q = < Py・Qz - Pz・Qy, Pz・Qx - Px・Qz, Px・Qy - Py・Qx >
		return Vector<N>{
			this->c[1] * rhs.c[2] - this->c[2] * rhs.c[1], 
			this->c[2] * rhs.c[0] - this->c[0] * rhs.c[2], 
			this->c[0] * rhs.c[1] - this->c[1] * rhs.c[0], 
		};
	}
 
	template <unsigned int N> 
	inline Vector<N - 1> Vector<N>::trunc() const {
		Vector<N - 1> res;
		copy(this->c, this->c + N - 1, res.c);
		return res;
	}
 
	template <unsigned int N> 
	inline ostream& operator <<(ostream& out, const Vector<N>& vec) {
		out << "{";
		for (auto c : vec.c) out << c << ", ";
		out << "}";
	}
};
 
#include <cmath>
 
namespace h3d {
	inline FLOAT ang2rad(const FLOAT& ang) {
		return M_PI * ang / 180.0;
	}
 
	inline bool nearEqual(const FLOAT& x, const FLOAT& y) {
		static const FLOAT THRESHOLD = 0.0000001;
		FLOAT dif = x - y;
		return dif > -THRESHOLD && dif < THRESHOLD;
	}
};
 
#include <iostream>
#include <string>
 
namespace h3d {
	using std::cout;
	using std::endl;
	using std::string;
 
	inline void assert(const string& pred_str, const bool& pred_res) throw(string) {
		if (pred_res) cout << "アサート成功: " << pred_str << endl;
		else throw "アサート失敗: " + pred_str;
	}
};
 
#include <memory>
 
namespace h3d {
	using std::shared_ptr;
 
	TestSet::TestSet() {
		this->tests.push_back(shared_ptr<Test>(new Test1));
		this->tests.push_back(shared_ptr<Test>(new Test2));
		this->tests.push_back(shared_ptr<Test>(new Test3));
	}
 
	void TestSet::run() const { for (auto iter : this->tests) iter->run(); }
};
 
#include <cmath>
 
namespace h3d {
	using std::sqrt;
 
	void Test1::run() const {
		Vector<3> p{2.0, 2.0, 1.0}, q{1.0, -2.0, 0.0};
		ASSERT(nearEqual(p.c[0], 2.0))
		ASSERT(nearEqual(p.c[1], 2.0))
		ASSERT(nearEqual(p.c[2], 1.0))
		ASSERT(nearEqual(q.c[0], 1.0))
		ASSERT(nearEqual(q.c[1], -2.0))
		ASSERT(nearEqual(q.c[2], 0.0))
 
		ASSERT(p.ext() == (Vector<4>{2.0, 2.0, 1.0, 1.0}))
		ASSERT(p.ext(0.0) == (Vector<4>{2.0, 2.0, 1.0, 0.0}))
		ASSERT(p.trunc() == (Vector<2>{2.0, 2.0}))
 
		ASSERT(p == (Vector<3>{2.0, 2.0, 1.0}))
		ASSERT(q == (Vector<3>{1.0, -2.0, 0.0}))
		ASSERT(p != (Vector<3>{1.0, -2.0, 0.0}))
		ASSERT(q != (Vector<3>{2.0, 2.0, 1.0}))
		ASSERT(-p == (Vector<3>{-2.0, -2.0, -1.0}))
		ASSERT(-q == (Vector<3>{-1.0, 2.0, 0.0}))
		ASSERT(p * 2.0 == (Vector<3>{4.0, 4.0, 2.0}))
		ASSERT(q * -3.0 == (Vector<3>{-3.0, 6.0, 0.0}))
		ASSERT(p / 2.0 == (Vector<3>{1.0, 1.0, 0.5}))
		ASSERT(q / -3.0 == (Vector<3>{1.0 / -3.0, -2.0 / -3.0, 0.0}))
 
		ASSERT(p + q == (Vector<3>{3.0, 0.0, 1.0}))
		ASSERT(q + p == (Vector<3>{3.0, 0.0, 1.0}))
		ASSERT(p - q == (Vector<3>{1.0, 4.0, 1.0}))
		ASSERT(q - p == (Vector<3>{-1.0, -4.0, -1.0}))
 
		Vector<3> r;
		r = p;
		ASSERT(r == (Vector<3>{2.0, 2.0, 1.0}))
		ASSERT((r *= 3.0) == (Vector<3>{6.0, 6.0, 3.0}))
		ASSERT((r /= -2.0) == (Vector<3>{-3.0, -3.0, -1.5}))
		ASSERT((r += (Vector<3>{5.0, -1.0, 0.0})) == (Vector<3>{2.0, -4.0, -1.5}))
		ASSERT((r -= (Vector<3>{0.0, -2.0, 4.0})) == (Vector<3>{2.0, -2.0, -5.5}))
		r.clear();
		ASSERT(nearEqual(r.c[0], 0.0))
		ASSERT(nearEqual(r.c[1], 0.0))
		ASSERT(nearEqual(r.c[2], 0.0))
 
		ASSERT(nearEqual(p.inpro(q), -2.0))
		ASSERT(p.outpro(q) == (Vector<3>{2.0, 1.0, -6.0}))
	}
};
 
#include <cmath>
 
namespace h3d {
	using std::cos;
	using std::sin;
	using std::sqrt;
 
	void Test2::run() const {
		ASSERT(nearEqual(Matrix<2>::IDENTITY.c[0][0], 1.0))
		ASSERT(nearEqual(Matrix<2>::IDENTITY.c[0][1], 0.0))
		ASSERT(nearEqual(Matrix<2>::IDENTITY.c[1][0], 0.0))
		ASSERT(nearEqual(Matrix<2>::IDENTITY.c[1][1], 1.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[0][0], 1.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[0][1], 0.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[0][2], 0.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[1][0], 0.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[1][1], 1.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[1][2], 0.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[2][0], 0.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[2][1], 0.0))
		ASSERT(nearEqual(Matrix<3>::IDENTITY.c[2][2], 1.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[0][0], 1.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[0][1], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[0][2], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[0][3], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[1][0], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[1][1], 1.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[1][2], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[1][3], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[2][0], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[2][1], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[2][2], 1.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[2][3], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[3][0], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[3][1], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[3][2], 0.0))
		ASSERT(nearEqual(Matrix<4>::IDENTITY.c[3][3], 1.0))
 
		Matrix<2> f{
			1.0,	-2.0, 
			-3.0,	4.0, 
		}, g{
			-9.0,	8.0, 
			7.0,	-6.0, 
		};
 
		ASSERT(nearEqual(f.c[0][0], 1.0))
		ASSERT(nearEqual(f.c[0][1], -2.0))
		ASSERT(nearEqual(f.c[1][0], -3.0))
		ASSERT(nearEqual(f.c[1][1], 4.0))
		ASSERT(nearEqual(g.c[0][0], -9.0))
		ASSERT(nearEqual(g.c[0][1], 8.0))
		ASSERT(nearEqual(g.c[1][0], 7.0))
		ASSERT(nearEqual(g.c[1][1], -6.0))
 
		ASSERT(f.ext() == (Matrix<3>{
			1.0,	-2.0,	0.0, 
			-3.0,	4.0,	0.0, 
			0.0,	0.0,	1.0, 
		}))
		ASSERT(f.ext((Matrix<3>{
			3.0,	-8.0,	5.0, 
			-2.0,	7.0,	4.0, 
			1.0,	-5.0,	0.0, 
		})) == (Matrix<3>{
			1.0,	-2.0,	5.0, 
			-3.0,	4.0,	4.0, 
			1.0,	-5.0,	0.0, 
		}))
		ASSERT(f.trunc() == (Matrix<1>{
			1.0, 
		}))
 
		ASSERT(f == (Matrix<2>{
			1.0,	-2.0, 
			-3.0,	4.0, 
		}))
		ASSERT(g == (Matrix<2>{
			-9.0,	8.0, 
			7.0,	-6.0, 
		}))
		ASSERT(f != (Matrix<2>{
			-9.0,	8.0, 
			7.0,	-6.0, 
		}))
		ASSERT(g != (Matrix<2>{
			1.0,	-2.0, 
			-3.0,	4.0, 
		}))
 
		ASSERT(-f == (Matrix<2>{
			-1.0,	2.0, 
			3.0,	-4.0, 
		}))
		ASSERT(-g == (Matrix<2>{
			9.0,	-8.0, 
			-7.0,	6.0, 
		}))
		ASSERT(f + g == (Matrix<2>{
			-8.0,	6.0, 
			4.0,	-2.0, 
		}))
		ASSERT(f - g == (Matrix<2>{
			10.0,	-10.0, 
			-10.0,	10.0, 
		}))
		ASSERT(f * g == (Matrix<2>{
			-23.0,	20.0, 
			55.0,	-48.0, 
		}))
		ASSERT(g * f == (Matrix<2>{
			-33.0,	50.0, 
			25.0,	-38.0, 
		}))
		ASSERT(f * Matrix<2>::IDENTITY == (Matrix<2>{
			1.0,	-2.0, 
			-3.0,	4.0, 
		}))
		ASSERT(Matrix<2>::IDENTITY * g == (Matrix<2>{
			-9.0,	8.0, 
			7.0,	-6.0, 
		}))
		ASSERT(f * 3.0 == (Matrix<2>{
			3.0,	-6.0, 
			-9.0,	12.0, 
		}))
		ASSERT(g / -2.0 == (Matrix<2>{
			4.5,	-4.0, 
			-3.5,	3.0, 
		}))
 
		Matrix<2> h;
		h = f;
		ASSERT(h == (Matrix<2>{
			1.0,	-2.0, 
			-3.0,	4.0, 
		}))
		ASSERT((h += (Matrix<2>{
			0.1,	-0.2, 
			-0.3,	0.4, 
		})) == (Matrix<2>{
			1.1,	-2.2, 
			-3.3,	4.4, 
		}))
		ASSERT((h -= (Matrix<2>{
			-0.9,	0.8, 
			0.7,	-0.6, 
		})) == (Matrix<2>{
			2.0,	-3.0, 
			-4.0,	5.0, 
		}))
		ASSERT((h *= (Matrix<2>{
			-4.0,	3.0, 
			-2.0,	1.0, 
		})) == (Matrix<2>{
			-2.0,	3.0, 
			6.0,	-7.0, 
		}))
		ASSERT((h *= 3.0) == (Matrix<2>{
			-6.0,	9.0, 
			18.0,	-21.0, 
		}))
		ASSERT((h /= -2.0) == (Matrix<2>{
			3.0,	-4.5, 
			-9.0,	10.5, 
		}))
		h.clear();
		ASSERT(nearEqual(h.c[0][0], 0.0))
		ASSERT(nearEqual(h.c[0][1], 0.0))
		ASSERT(nearEqual(h.c[1][0], 0.0))
		ASSERT(nearEqual(h.c[1][1], 0.0))
 
		ASSERT(nearEqual((Matrix<2>{
			2.0,	7.0, 
			-3.0,	0.5, 
		}).det(), 22.0))
		ASSERT(nearEqual((Matrix<3>{
			0.0, 0.0, 1.0, 
			0.0, 1.0, 0.0, 
			1.0, 0.0, 0.0, 
		}).det(), -1.0))
		ASSERT(nearEqual((Matrix<3>{
			0.5,				sqrt(3.0) / 2.0,	0.0, 
			-sqrt(3.0) / 2.0,	0.5,				0.0, 
			0.0,				0.0,				1.0, 
		}).det(), 1.0))
		ASSERT(nearEqual((Matrix<3>{
			5,	7,	1, 
			17,	2,	64, 
			10,	14,	2, 
		}).det(), 0.0))
 
		ASSERT((Matrix<3>{
			2.0,	0.0,	0.0, 
			0.0,	3.0,	0.0, 
			0.0,	0.0,	4.0, 
		}).inv() == (Matrix<3>{
			12.0,	0.0,	0.0, 
			0.0,	8.0,	0.0, 
			0.0,	0.0,	6.0, 
		}) / 24.0)
		ASSERT((Matrix<3>{
			2.0,	0.0,	0.0, 
			0.0,	3.0,	0.0, 
			0.0,	0.0,	4.0, 
		}) * ((Matrix<3>{
			12.0,	0.0,	0.0, 
			0.0,	8.0,	0.0, 
			0.0,	0.0,	6.0, 
		}) / 24.0) == Matrix<3>::IDENTITY)
		ASSERT((Matrix<3>{
			1.0,	0.0,	0.0, 
			0.0,	2.0,	2.0, 
			3.0,	0.0,	8.0, 
		}).inv() == (Matrix<3>{
			16.0,	0.0,	0.0, 
			6.0,	8.0,	-2.0, 
			-6.0,	0.0,	2.0, 
		}) / 16.0)
		ASSERT((Matrix<3>{
			1.0,	0.0,	0.0, 
			0.0,	2.0,	2.0, 
			3.0,	0.0,	8.0, 
		}) * ((Matrix<3>{
			16.0,	0.0,	0.0, 
			6.0,	8.0,	-2.0, 
			-6.0,	0.0,	2.0, 
		}) / 16.0) == Matrix<3>::IDENTITY)
		FLOAT theta_rad = ang2rad(60.0);
		FLOAT cos_res = cos(theta_rad);
		FLOAT sin_res = sin(theta_rad);
		ASSERT((Matrix<3>{
			cos_res,	0.0,	-sin_res, 
			0.0,	1.0,	0.0, 
			sin_res,	0.0,	cos_res, 
		}).inv() == (Matrix<3>{
			cos_res,	0.0,	sin_res, 
			0.0,	1.0,	0.0, 
			-sin_res,	0.0,	cos_res, 
		}))
		ASSERT((Matrix<3>{
			cos_res,	0.0,	-sin_res, 
			0.0,	1.0,	0.0, 
			sin_res,	0.0,	cos_res, 
		}) * (Matrix<3>{
			cos_res,	0.0,	sin_res, 
			0.0,	1.0,	0.0, 
			-sin_res,	0.0,	cos_res, 
		}) == Matrix<3>::IDENTITY)
		ASSERT((Matrix<4>{
			1.0,	0.0,	0.0,	4.0, 
			0.0,	1.0,	0.0,	3.0, 
			0.0,	0.0,	1.0,	7.0, 
			0.0,	0.0,	0.0,	1.0, 
		}).inv() == (Matrix<4>{
			1.0,	0.0,	0.0,	-4.0, 
			0.0,	1.0,	0.0,	-3.0, 
			0.0,	0.0,	1.0,	-7.0, 
			0.0,	0.0,	0.0,	1.0, 
		}))
		ASSERT((Matrix<4>{
			1.0,	0.0,	0.0,	4.0, 
			0.0,	1.0,	0.0,	3.0, 
			0.0,	0.0,	1.0,	7.0, 
			0.0,	0.0,	0.0,	1.0, 
		}) * (Matrix<4>{
			1.0,	0.0,	0.0,	-4.0, 
			0.0,	1.0,	0.0,	-3.0, 
			0.0,	0.0,	1.0,	-7.0, 
			0.0,	0.0,	0.0,	1.0, 
		}) == Matrix<4>::IDENTITY)
 
		ASSERT((Matrix<2>{
			1.0, -3.0, 
			-4.0, 6.0, 
		}) * (Vector<2>{2.0, 5.0}) == (Vector<2>{-13.0, 22.0}))
		ASSERT((Matrix<2>{
			0.0, 8.0, 
			-9.0, 3.0, 
		}) * (Vector<2>{1.0, -4.0}) == (Vector<2>{-32.0, -21.0}))
	}
};
 
#include <cmath>
 
namespace h3d {
	using std::sqrt;
 
	void Test3::run() const {
		ASSERT(scaling((Vector<3>{1.0, -2.0, 0.0})) == (Matrix<3>{
			1.0,	0.0,	0.0, 
			0.0,	-2.0,	0.0, 
			0.0,	0.0,	0.0, 
		}))
		ASSERT(scaling((Vector<3>{0.0, 9.0, -8.0})) == (Matrix<3>{
			0.0,	0.0,	0.0, 
			0.0,	9.0,	0.0, 
			0.0,	0.0,	-8.0, 
		}))
 
		ASSERT(rotation((Vector<3>{1.0, 0.0, 0.0}), ang2rad(30.0)) == (Matrix<3>{
			1.0,				0.0,				0.0, 
			0.0,				sqrt(3.0) / 2.0,	-0.5, 
			0.0,				0.5,				sqrt(3.0) / 2.0, 
		}))
		ASSERT(rotation((Vector<3>{0.0, 1.0, 0.0}), ang2rad(-45.0)) == (Matrix<3>{
			1.0 / sqrt(2.0),	0.0,				-1.0 / sqrt(2.0), 
			0.0,				1.0,				0.0, 
			1.0 / sqrt(2.0),	0.0,				1.0 / sqrt(2.0), 
		}))
		ASSERT(rotation((Vector<3>{0.0, 0.0, 1.0}), ang2rad(60.0)) == (Matrix<3>{
			0.5,				-sqrt(3.0) / 2.0,	0.0, 
			sqrt(3.0) / 2.0,	0.5,				0.0, 
			0.0,				0.0,				1.0, 
		}))
 
		ASSERT(translation((Vector<3>{1.0, -2.0, 0.0})) == (Matrix<4>{
			1.0,	0.0,	0.0,	1.0, 
			0.0,	1.0,	0.0,	-2.0, 
			0.0,	0.0,	1.0,	0.0, 
			0.0,	0.0,	0.0,	1.0, 
		}))
		ASSERT(translation((Vector<3>{0.0, 9.0, -8.0})) == (Matrix<4>{
			1.0,	0.0,	0.0,	0.0, 
			0.0,	1.0,	0.0,	9.0, 
			0.0,	0.0,	1.0,	-8.0, 
			0.0,	0.0,	0.0,	1.0, 
		}))
	}
};
 
#include <iostream>
#include <string>
 
int main() {
	using std::cerr;
	using std::endl;
	using std::string;
 
	try {
		h3d::TestSet().run();
	}
	catch (const string& msg) {
		cerr << msg << endl;
		return 1;
	}
	return 0;
}