「prolog勉強プロジェクトオイラー91~100」の編集履歴(バックアップ)一覧に戻る
prolog勉強プロジェクトオイラー91~100」を以下のとおり復元します。
プロジェクトオイラーの問題を堀江伸一さんがPrologで解くページ。
*問い91
原点が特殊処理が必要なくらいで後は点OAのベクトルを90度回して最小公倍数でその長さ割ったベクトルの定数倍になる点を数えればOK。

 gcd(0, B, B).
 gcd(A, B, G) :- A > 0, R is B mod A, gcd(R, A, G). 
 
 min(X,Y,X):-X<Y,!.
 min(_,Y,Y):-!.
 
 point_count(_,0,1000):-!.
 point_count(L,D,Result):-
 	Result is L//D.
 
 r90(0,_,1,0):-!.
 r90(_,0,0,1):-!.
 r90(DY,DX,DY1,DX1):-
 	gcd(DY,DX,G),
 	DY1 is DX//G,
 	DX1 is DY//G.
  
 calc_perm(Y,X,Ans,Ans1):-
 	r90(Y,X,DY,DX),
 	X1 is 50-X,
 	Y1 is 50-Y,
	point_count(X ,DX,T1),
 	point_count(Y1,DY,T2),
  	point_count(X1,DX,T3),
 	point_count(Y,DY,T4),
 	min(T1,T2,Re1),
 	min(T3,T4,Re2),
 	Ans1 is Ans+Re1+Re2.
 
 search(51,0,Ans):-!,write(Ans).
 search(Y,51,Ans):-
 	!,
 	Y1 is Y+1,
 	search(Y1,0,Ans).
 search(Y,X,Ans):-
 	calc_perm(Y,X,Ans,Ans1),
 	X1 is X+1,
 	search(Y,X1,Ans1).
 
 main91:-search(0,1,2500).




*問い92 Square digit chains
http://projecteuler.net/problem=92
この問題はC++なら配列を使えば結構な速度で答えが出る。
Prologでは遅い。
少しだけ工夫が必要だ問題の関数をfと定義すると。
f(102345)=f(201345)
のように数字を並べ替えたものは同じ値になることとどうやら
f(1~9999999)<=567
となることを利用して計算量を下げれそうに見える。
これで試してみることにしよう。
計算を簡単にするために0012407のような0づめを許して組み合わせ数を求めることにする。
動作が遅いPrologだがほとんど瞬きする間に答えが出た。
こんなものだろう。
C++に翻訳して試しに一千万を1兆に変えてみたところやはりまたたきする間に答えが出た。


 to_perm(N,Result):-nth0(N,[1,1,2,6,24,120,720,5040,40320],Result).
 search(Num,7,Div,NowDiv,_,[Num,Result]):-
 	!,
 	Num>0,
 	to_perm(NowDiv,NowDiv1), 
 	Result is 5040//(Div*NowDiv1).
 search(Num,Deep,Div,NowDiv,[N|Nums],Result):-
 	NowDiv1 is NowDiv+1,
 	Num1 is Num+N*N,
 	Deep1 is Deep+1,
 	search(Num1,Deep1,Div,NowDiv1,[N|Nums],Result).
 
 search(Num,Deep,Div,NowDiv,[_|Nums],Result):-
 	to_perm(NowDiv,NowDiv1),
 	Div1 is Div*NowDiv1,
 	search(Num,Deep,Div1,0,Nums,Result).
 
 sum([],[N,P],[[N,P]]):-!.
 sum([[N,P]|Rest],[N,P1],Result):-
 	!,
 	P2 is P+P1,
 	sum(Rest,[N,P2],Result).
 sum([[N,P]|Rest],[N1,P1],[[N1,P1]|Result]):-
 	sum(Rest,[N,P],Result).
 
 calc(0,Sum,Sum):-!.
 calc(N,Sum,Result):-
 	M is N mod 10,
 	Sum1 is Sum+M*M,
 	N1 is N//10,
 	!,
 	calc(N1,Sum1,Result).
 
 search(1):-!,fail.
 search(89):-!.
 search(N):-
 	calc(N,0,N1),
 	search(N1).
 
 seed(N):-
 	between(1,567,N),
 	search(N).
 
 sum2([],_,Sum,Sum):-!.
 sum2(_,[],Sum,Sum):-!.
 sum2([[N,P]|Perm],[N|OKs],Sum,Result):-
 	!,
 	Sum1 is Sum+P,
 	sum2(Perm,OKs,Sum1,Result).
 sum2([[N,_]|Perm],[N1|OKs],Sum,Result):-
 	N<N1,
 	!,
 	sum2(Perm,[N1|OKs],Sum,Result).
 sum2(Perm,[_|OKs],Sum,Result):-
 	!,
 	sum2(Perm,OKs,Sum,Result).
 
 main92:-Nums=[0,1,2,3,4,5,6,7,8,9],
 	findall(P,search(0,0,1,0,Nums,P),Perm),
 	msort(Perm,Perm1),
 	[[N,_]|_]=Perm1,
 	sum(Perm1,[N,0],Perm2),
 	findall(N1,seed(N1),OKs),
 	sum2(Perm2,OKs,0,Ans),write(Ans).





*問い92 Bcc5.5 C++ Code
1000万を一兆にしてみて計算、またたきする間に答えが出る。
答えは 850878696414通り


 #include<stdio.h>
 #include<iostream>
 #include<vector> 
 
 std::vector<__int64> to_perm;
 int oks[1001]={0};
 const int up=1000;
 const int keta=12;//7 1000万 8 1億 9 10億 10 100億 11 10000億、12 一兆
 
 int calc(int n){
 	if(n==0){
 		return 0;
 	}else{
 		int m=n % 10;
 		return calc(n/10)+m*m;
 	}
 } 
  
 int calc_chain(int n){
 	if(oks[n]==0){
 		oks[n]=calc_chain(calc(n));
 	}
 	return oks[n];
 } 
  
 __int64 search(int num,int deep,__int64 div,int nowDiv,int p){
 	__int64 ans=0;
 	if(deep==keta){
 		return to_perm[deep]/(div*to_perm[nowDiv])*oks[num];
 	}else{
 		ans+=search(num+p*p,deep+1,div,nowDiv+1,p);
 		if(p<9){
 			ans+=search(num,deep,div*to_perm[nowDiv],0,p+1);
 		}
 	}
 	return ans;
 }
 
 int main(){
 	oks[89]=1;
 	oks[1]=-1;
  	for(int i=1;i<=up;i++){
 		calc_chain(i);
 	}
 	for(int i=1;i<=up;i++){
 		if(oks[i]!=1)oks[i]=0;
 	}
 	__int64 p=1;
  	to_perm.push_back(p);
 	for(int i=1;i<=12;i++){
 		p*=i;
 		to_perm.push_back(p);
 	}
 	__int64 ans=search(0,0,1,0,0);
 	std::cout<<ans;
 }



*Problem 94 「擬正三角形」 †
http://odz.sakura.ne.jp/projecteuler//index.php?cmd=read&page=Problem%2094
色々ショートカットテクニックを考えたのですが結局原始ピタゴラス数でNを1ずつ増加させて、Nから斜辺と底辺*2が1差になるMを2次方程式の解法から求めて解きました。
計算の終了条件は少し手抜きです。
 calc(N,Len):-述語の適用式のmember化は少し楽しいテクニックを思いつけたと思っています。
多分Prologをわかってる人からみたら初歩だと思いますが自力でこういうのを思いつけたとこの楽しさは格別ですね。





以下コード。
 over(N):-1000000000<N.
 is_square(N,Result):-sqrt(N,N1),
 	Result is floor(N1), 
 	N=:=Result*Result.
 gcd(0, B, B).
 gcd(A, B, G) :- A > 0, R is B mod A, gcd(R, A, G).
 
 add([],Ans,Ans):-!.
 add([X|Lens],Ans,Result):-
 	!,
 	Ans1 is Ans+X,
 	add(Lens,Ans1,Result).
 
 calc(N,Len):-
 	member([S1,S2,S3]
 		     ,[[3*N*N-1,2*N+T1,2*M*M+2*N*N+4*M*N],
 		       [3*N*N+1,2*N+T1,2*M*M+2*N*N+4*M*N],
  		       [3*N*N+1,T1,4*M*M],
 		       [3*N*N-1,T1,4*M*M]]),
 	T is S1,
 	is_square(T,T1),
 	M is S2,
 	M>N,
 	(M-N) mod 2=:=1,
 	gcd(M,N,G),
 	G=:=1,
 	Len is S3,
 	not(over(Len)).
 
 search(N,Ans,Ans):-
 	T is 4*N*N,over(T),!.
 search(N,Ans,Result):-
 	N1 is N+1,
 	findall(Len,calc(N,Len),Lens),
 	add(Lens,Ans,Ans1),
 	search(N1,Ans1,Result).
 main94:-search(1,0,Ans),write(Ans).




*問い95
この問題はC++なら中学生でも解ける問題である。
Prologで高速に解くとなるとかなりの工夫が必要に見える。
正直言うと解き方を全く思いつかない。
原理的に言って無いのかもしれない。

復元してよろしいですか?