プロジェクトオイラー問177

http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%20177
整数角の四角形の個数を数える問題。

詳細はリンク先を参照してください。
ややこしい問題なので下記の内容を紙に図に書いて考えてください。

角が全て整数の三角形ABCの集合Sを考える。
ABの長さを1とするとSは有限となる。
Sから一つずつ選びE1とする。
すると個々の三角形において正弦定理よりACの長さが出る。
角が全て整数の三角形ACDの集合Uを考える。
Uのうち一つずつ選びE2とする。

E1とE2を張り合わせた四角形ABCDとしこれをE3と考える。
ACの長さがでているのでE1と張り合わせられるE2はACの長さがE1のACの長さと同じでなくてはいけない。
これによりUは個数が有限となる。
ACの長さからADとDCの長さが正弦定理より算出される。

整数角の三角形を張り合わせているので8つの角のうち4つは整数角であることが保障される。

余弦定理よりDBの長さが算出できる。

AB,BD、DBの長さが判明するので余弦定理より
Cos角ADB
Cos角ABDと5個めと6個めの角度が算出できる。
これらが整数角のCosとほぼ同じ値かつ角ABC、角ADCより小さい値ならE3は条件を満たす四角形である。
なぜなら残る角度は整数-整数だから整数である。

後は計算できた8つの角を使って四角形の重複数え上げを排除して数えれば答えとなります。
計算に数分かかりますが簡単な枝狩りを幾つかすればこのコードでも数倍くらいは速くなるはずです。

実装と重複数え上げの排除で細かいミスをして試行錯誤の時間を含めてコードを書くのに3時間くらいかかりました。

#include<stdio.h>
#include<math.h>
#include<map>
#include<algorithm>
#include<vector>
#include<set>


std::set<double> memo;
double LIMIT=0.000000001;


struct E2{
	int rs[8];
	
	void set(const int as[8]){
		E2 e2,e2A;
		e2.set2(as);
		e2A.set2(as);
 		
		for(int i=0;i<4;i++){
			if(e2A.hikaku(e2)==true)e2A.set2(e2.rs);
			e2.round();
		}
		e2.reverse();
		for(int i=0;i<4;i++){
			if(e2A.hikaku(e2)==true)e2A.set2(e2.rs);
			e2.round();
		}
		set2(e2A.rs);
	}
	void set2(const int as[8]){
		for(int i=0;i<8;i++)rs[i]=as[i];
	}
	bool hikaku(const E2& e)const{
		for(int i=0;i<8;i++){
			if(rs[i]!=e.rs[i])return rs[i]<e.rs[i];
		}
		return false;
	}
	
	void round(){
		int r6=rs[6];
		int r7=rs[7];
		for(int i=7;i>=2;i--){
			rs[i]=rs[i-2];
		}
		rs[1]=r7;
		rs[0]=r6;
	}
	void reverse(){
		std::swap(rs[0],rs[1]);
		std::swap(rs[7],rs[2]);
		std::swap(rs[6],rs[3]);
		std::swap(rs[4],rs[5]);
	}
	
};

struct E{
	int rs[8];
	bool ok(){
		for(int i=0;i<8;i++){
			if(rs[i]<1)return false;
		}
		for(int i=0;i<8;i+=2){
			if(rs[i]+rs[i+1]>=180)return false;
		}
		
		return true;
	}	
	void set(int cad,int cab,int dba,int dbc,int acb,int acd,int cdb,int bda){
		rs[0]=cad;
		rs[1]=cab;
		rs[2]=dba;
		rs[3]=dbc;
		rs[4]=acb;
		rs[5]=acd;
		rs[6]=cdb;
		rs[7]=bda;
		E2 e2;
		e2.set(rs);
		for(int i=0;i<8;i++)rs[i]=e2.rs[i];
	}
	bool operator<(const E& e1)const{
		for(int i=0;i<8;i++){
			if(rs[i]!=e1.rs[i])return rs[i]<e1.rs[i];
		}
		return false;
	}
};



struct D{
	int a,b,c;
};
int kakudo(double r){
	//-1なら外れ、
	std::set<double>::iterator it;
	if(r<=-1||1<=r)return -1;
	r=acos(r)/M_PI*180.0;
	it=memo.lower_bound(r);
	double r2;
	if(it!=memo.end()){
		r2=(*it);
		if(fabs(r-r2)<=LIMIT)return (int)r2;
	}
	if(it!=memo.begin()){
		it--;
		r2=(*it);
		if(fabs(r-r2)<=LIMIT)return (int)r2;
	}
	return -1;
} 

double seigen(double c,double rC,double rB){
	double sinC=sin(rC/180.0*M_PI);
	double sinB=sin(rB/180.0*M_PI);
	return c*sinB/sinC;
}
double yogen(double b,double c,double rA){
	return sqrt(b*b+c*c-2*b*c*cos(rA/180.0*M_PI));
}
double calcCos(double a,double b,double c){
	return (b*b+c*c-a*a)/(2*b*c);
} 

 
int main(){
	for(int i=1;i<180;i++){
		memo.insert(i);
	}
	std::set<E> count;
	std::vector<D> vec;
	D d1,d2;
	for(int a=1;a<180;a++){
		for(int b=1;a+b<180;b++){
			int c=180-a-b;
			d1.a=a;
			d1.b=b;
			d1.c=c;
			vec.push_back(d1);
		}
	}
	for(int i=0;i<vec.size();i++){
		for(int j=0;j<vec.size();j++){
			d1=vec[i];
 			d2=vec[j];
			if(d1.b+d2.b>=180)continue;
			if(d1.c+d2.c>=180)continue;
			
			double lenAC=seigen(1,d1.c,d1.a);
			double lenAD=seigen(lenAC,d2.a,d2.c);
			double lenBD=yogen(lenAD,1,d1.b+d2.b);
			double cosADB=calcCos(1,lenBD,lenAD);
			double cosABD=calcCos(lenAD,1,lenBD);
			int rADB=kakudo(cosADB);
			int rABD=kakudo(cosABD);
			
			int rDAC=d2.b;
			int rCAB=d1.b;
			rABD=rABD;
			int rDBC=d1.a-rABD;
 			int rACB=d1.c;
			int rACD=d2.c;
			int rBDC=d2.a-rADB;
			rADB=rADB;
			E e1;
 			e1.set(rDAC,rCAB,rABD,rDBC,rACB,rACD,rBDC,rADB);
			if(e1.ok()==false)continue;
		
 			
			count.insert(e1);
		}
 	}
	printf("ans=%d\n",count.size());
}
最終更新:2015年11月02日 17:09