「prolog勉強プロジェクトオイラー21~30」の編集履歴(バックアップ)一覧に戻る

prolog勉強プロジェクトオイラー21~30 - (2014/01/07 (火) 02:07:53) のソース

プロジェクトオイラーの問題をPrologで解くページ。
創価学会の皆さまから小学校の算数までしかできないと言われている堀江伸一さんがこのページのコードを書いているようです。




*Problem 21 「階乗の数字和」 †
http://projecteuler.net/problem=21
10000以下の友愛数の和を求めよという問題。
そのまま定義通り計算するだけです。
複数の関数で一つの関数を表現するPrologて独特。

 yakusuu_sum(1,_,1,Sum,Sum):-!.
 yakusuu_sum(N,M,1,Sum,Sum1):-N<M*M,!,Sum1 is Sum*(N+1).
 yakusuu_sum(N,M,Multi,Sum,Result):-
	0=:=N mod M,!,
	N1 is N//M,
	Multi1 is Multi*M+1,
	yakusuu_sum(N1,M,Multi1,Sum,Result).
 yakusuu_sum(N,M,Multi,Sum,Result):-
	Sum1 is Sum*Multi,
	M1 is M+1,
	yakusuu_sum(N,M1,1,Sum1,Result).
 yakusuu_sum_w(N,Result):- 
	yakusuu_sum(N,2,1,1,Result).
 
 check(N):-yakusuu_sum_w(N,N1),N2 is N1-N,N\==N2,yakusuu_sum_w(N2,N3),
         N4 is N3-N2,N=:=N4.
 calc(N):-
 	between(2,10000,N),check(N).
 sum([],Sum,Sum):-!.
 sum([X|Rest],Sum,Result):-Sum1 is Sum+X,sum(Rest,Sum1,Result).
 
 main21:-findall(N,calc(N),List),sum(List,0,Ans),write(Ans).




*Problem 22 「名前のスコア」 †
http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2022
人名の書かれたテキストファイルを読み込み、スコアに変換してその合計をこたえよという問題。

コード
Prologの徹底ぶりを思い知った問題。
まさかファイル読み込み関数までバックトラックの対象だとは思わず、他の言語と同じノリで一回読んだらそれで終わり的な発想でコードを書いてすこしだけはまった。
気がついて慌てて修正。
このコードは正しいのだろうか?
一応末尾再帰で再帰の深さを回避はしてるつもり。
findallを使う方が正しくないだろうか?


 last(-1).
 spliter(34).
 
 read_name(Names,Name):-
	get0(C),
	(spliter(C)->
        reverse(Name,Name1),karayomi_read([Name1|Names]);
	C1 is C-"A"+1,
	read_name(Names,[C1|Name])).
 
 karayomi_read(Names):-get0(C),karayomi(Names,C).
 karayomi(Names,C):-last(C),!,calc(Names).
 karayomi(Names,C):-spliter(C),!,read_name(Names,[]).
 karayomi(Names,_):-karayomi_read(Names).
 
 
 sum_score([],Sum,Sum).
 sum_score([X|Rest],Sum,Result):-Sum1 is Sum+X,sum_score(Rest,Sum1,Result).
 
 score_calc([],_,Score,Score).
 score_calc([Name|Rest],N,Score,Result):-sum_score(Name,0,S),Score1 is Score+N*S,
	N1 is N+1,
	score_calc(Rest,N1,Score1,Result).
 
 
 calc(Names):-seen,sort(Names,Names1),
	score_calc(Names1,1,0,Ans),write(Ans).
 
 main22:-seen,see('e22.txt'),karayomi_read([]).


*問い23 Non-abundant sums
http://projecteuler.net/problem=23
完全数とは, その数の真の約数の和がそれ自身と一致する数のことである. たとえば, 28の真の約数の和は, 1 + 2 + 4 + 7 + 14 = 28 であるので, 28は完全数である.
真の約数の和がその数よりも少ないものを不足数といい, 真の約数の和がその数よりも大きいものを過剰数と呼ぶ.
12は, 1 + 2 + 3 + 4 + 6 = 16 となるので, 最小の過剰数である. よって2つの過剰数の和で書ける最少の数は24である. 数学的な解析により, 28123より大きい任意の整数は2つの過剰数の和で書けることが知られている. 2つの過剰数の和で表せない最大の数がこの上うよりも小さいことは分かっているのだが, この上限を減らすことが出来ていない.
2つの過剰数の和で書き表せない正の整数の総和を求めよ.


コード
配列のある言語だったら簡単で高速に解ける問題だけど、Prologには配列がないためにちょっとしたパズルになってしまう問題。
一応末尾再帰で再帰にまつわる問題を回避しているので、再帰が深くなりすぎて停止することはないはず。
末尾に置いた述語はバックトラックが起きないよう述語内でカットしておかないと末尾再帰にならないようである。
このへんPrologの利点が欠点になってる。
末尾再帰になるよう色々プログラマの方でチェックしてあげないといけない。


 yakusuu_sum(N,M,Multi,Sum,Result):-
 	0=:=N mod M,!,
 	N1 is N//M,
 	Multi1 is Multi*M+1,
 	yakusuu_sum(N1,M,Multi1,Sum,Result).
 yakusuu_sum(N,M,Multi,Sum,Result):-
 	!,Sum1 is Sum*Multi,
 	M1 is M+1,
 	yakusuu_sum(N,M1,1,Sum1,Result).
 yakusuu_sum_w(N,Result):-
 	yakusuu_sum(N,2,1,1,Result).
 sa_list(_,[],[]):-!.
 sa_list(N,[X|Rest],[N1|Result]):-
 	N1 is N-X,
 	sa_list(N,Rest,Result).
 check([],_,List,N,Ans):-!,add_next(List,N,Ans,N).
 check(_,[],List,N,Ans):-!,add_next(List,N,Ans,N).
 check([X1|_],[X1|_],List,N,Ans):-!,add_next(List,N,Ans,0).
 check([X1|Rest1],[X2|Rest2],List,N,Ans):-X1<X2,!,
 	check(Rest1,[X2|Rest2],List,N,Ans).
 check([X1|Rest1],[X2|Rest2],List,N,Ans):-X1>X2,
 	check([X1|Rest1],Rest2,List,N,Ans).
 check_w(List,List2,N,Ans):-
 	reverse(List,List1),!,
 	check(List1,List2,List,N,Ans). 
 add_next(List,N,Ans,Add):-
 	yakusuu_sum_w(N,Re),
 	N1 is N+1,
 	N2 is N*2,
 	Re>N2,!,
 	Ans1 is Ans + Add,
 	get_kazyou_list([N|List],N1,Ans1).
 add_next(List,N,Ans,Add):-
 	N1 is N+1,!,
 	Ans1 is Ans+Add,
 	get_kazyou_list(List,N1,Ans1).
 get_kazyou_list(List,28123,Ans):-!,nl,write([Ans]),length(List,Len),write(Len).
 get_kazyou_list(List,N,Ans):-!,
 	sa_list(N,List,List1),
 	check_w(List,List1,N,Ans).
 main23:-get_kazyou_list([],1,0).






*Problem 24 「辞書式順列」 †
http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2024
順列とはモノの順番付きの並びのことである. たとえば, 3124は数 1, 2, 3, 4 の一つの順列である. すべての順列を数の大小でまたは辞書式に並べたものを辞書順と呼ぶ. 0と1と2の順列を辞書順に並べると
012 021 102 120 201 210
になる.
0,1,2,3,4,5,6,7,8,9からなる順列を辞書式に並べたときの100万番目はいくつか?


手計算で解いた方が100倍早い問題。
一行で書ける言語で解くのが正しいが、これをProlog的に書くとこんな風になる。
この問題をプログラムで解くこと自体がナンセンスと言われても仕方ない。


 facts(0,[1]):-!.
 facts(X,[T1|Result]):-
 	X1 is X-1,
 	facts(X1,Result),
 	[T|_]=Result,
 	T1 is T*X.
 calc(0,[],_,Ans):-!,reverse(Ans,Ans1),write(Ans1).
 calc(X,[M|Facters],Nums,Ans):-
 	X1 is X mod M,
 	P1 is X //M,
 	nth0(P1,Nums,N1),
 	select(N1,Nums,Nums1),
 	calc(X1,Facters,Nums1,[N1|Ans]). 
 main24:-
 	facts(9,Facters),write(Facters),calc(999999,Facters,[0,1,2,3,4,5,6,7,8,9],[]).


*問い25 1000-digit Fibonacci number
http://projecteuler.net/problem=25
Problem 25 「1000桁のフィボナッチ数」 †
フィボナッチ数列は以下の漸化式で定義される:

Fn = Fn-1 + Fn-2, ただし F1 = 1, F2 = 1.
最初の12項は以下である.

F1 = 1
F2 = 1
F3 = 2
F4 = 3
F5 = 5
F6 = 8
F7 = 13
F8 = 21
F9 = 34
F10 = 55
F11 = 89
F12 = 144
12番目の項, F12が3桁になる最初の項である.
1000桁になる最初の項の番号を答えよ.
解法
多桁計算を標準でサポートしている言語ではこれは簡単な問題。
まあフィナボッチ数列の公式を使って桁数を求めるのも手です。
賢い人はそっちで計算してください、私は賢くないので単純に計算しました。

 fibo(A,_,N):-10^999=<A,!,write(N).
 fibo(A,B,N):-C is A+B,N1 is N+1,fibo(B,C,N1).
 main25:-fibo(1,1,1).






*問い26 Reciprocal cycles
http://projecteuler.net/problem=26

単位分数とは分子が1の分数である. 分母が2から10の単位分数を10進数で表記すると次のようになる.
1/2 = 0.5
1/3 = 0.(3)
1/4 = 0.25
1/5 = 0.2
1/6 = 0.1(6)
1/7 = 0.(142857)
1/8 = 0.125
1/9 = 0.(1)
1/10 = 0.1
0.1(6)は 0.166666... という数字であり, 1桁の循環節を持つ. 1/7 の循環節は6桁ある.
d < 1000 なる 1/d の中で小数部の循環節が最も長くなるような d を求めよ.
解法

一番素朴な方法で実装。
もうちょっとましな処理を書きたいところだなこれは。
この手法では小学生と変わらない。
そのうえ割り算のたびにリストを精査してるから遅い。
せめて配列があれば楽なんだけど。
9、、、9で割ればよいのですが循環小数の循環節開始点の判別が出来ないのが問題。

http://www2r.biglobe.ne.jp/kosanhp/math/junkan.pdf
この資料が役立つらしい。


 zyunn_len(_,Num,_,M,0):-0=:=Num mod M,!.
 zyunn_len(List,Num,P,M,Ans):-
 	Num<M,!,Num1 is Num*10,zyunn_len(List,Num1,P,M,Ans).
 zyunn_len(List,Num,P,M,Ans):-
 	P1 is P+1,
 	Num1 is Num mod M,
 	(member([X,Num1],List)->
 	Ans is P-X;
 	zyunn_len([[P,Num1]|List],Num1,P1,M,Ans)).
 
 calc(Num,Len):-
 	between(1,1000,Num),
 	zyunn_len([],1,0,Num,Len).
 
 max([],[0,0]):-!.
 max([[Num,Len]|List],[Num2,Len2]):-
 	max(List,[Num1,Len1]),
 	(Len1<Len->
 	Len2 is Len ,Num2 is Num;
 	Len2 is Len1,Num2 is Num1).
 
 main26:-findall([Num,Len],calc(Num,Len),List),
 	max(List,Ans),write(Ans).



他人のCコードを参考に書きなおしたコード。
循環節がどこから始まるにせよ小数点以下N桁目よりは開始点が手前にあるという性質を利用して、N桁目まで求め、後はその時の余りが出てくるまでさらに割り続ければ循環節の長さが出てくる。
コードを見てるとなるほどと思う発想。
こういうちょっとした発想、思いつけるようになったら仕事でもちょっとは使えそうな気がする。

 calc(0,_,_,_,0):-!.
 calc(N,M,R,Len,Len1):-N1 is (N*10) mod M,N1=:=R,!,Len1 is Len+1.
 calc(N,M,R,Len,Result):-
 	Len1 is Len+1,
 	N1 is (N*10) mod M,
 	calc(N1,M,R,Len1,Result).
 
 calc_mod(N,R,N,R):-!.
 calc_mod(N,R,C,Result):-
 	C1 is C+1,
 	R1 is (R*10) mod N,
 	calc_mod(N,R1,C1,Result).
 
 search([Len,N]):-
 	between(2,1000,N),
 	calc_mod(N,1,0,R),
 	calc(R,N,R,0,Len).
 main26:-findall([Len,N],search([Len,N]),List),
 	sort(List,List1),write(List1).



*問い27 Quadratic primes
http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2027
オイラーは以下の二次式を考案している:
n^2 + n + 41.
この式は, n を0から39までの連続する整数としたときに40個の素数を生成する. しかし, n = 40 のとき 402 + 40 + 41 = 40(40 + 1) + 41 となり41で割り切れる. また, n = 41 のときは 41^2 + 41 + 41 であり明らかに41で割り切れる.
計算機を用いて, 二次式 n^2 - 79n + 1601 という式が発見できた. これは n = 0 から 79 の連続する整数で80個の素数を生成する. 係数の積は, -79 × 1601 で -126479である.
さて, |a| < 1000, |b| < 1000 として以下の二次式を考える (ここで |a| は絶対値): 例えば |11| = 11 |-4| = 4である。
n^2 + an + b
n = 0 から始めて連続する整数で素数を生成したときに最長の長さとなる上の二次式の, 係数 a, b の積を答えよ.

解法
色々考えたが、処理手順を並べ替えるだけの素朴な発想では計算手数が落ちない。
高度な数学知識を使った方法は知らないので結論として素朴に全部試している。


 not_prime(N):-N<2,!,true.
 not_prime(N):-
 	sqrt(N,A),
  	B is floor(A),
 	between(2,B,M),
 	M1 is N mod M,
 	0=:=M1,!,true.
 
 calc(N,A,B,N):-V is N*N+A*N+B,not_prime(V),!.
 calc(N,A,B,Result):-N1 is N+1,calc(N1,A,B,Result).
 
 search(N,A,B):-
 	between(2,1000,B),
 	not(not_prime(B)),
 	between(-1000,1000,A),
 	calc(0,A,B,N).
 
 main27:-findall([N,A,B],search(N,A,B),List),
 	sort(List,List1),write(List1).




*Problem 28 「螺旋状に並んだ数の対角線」 †
http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2028
格子の中でらせん状に生成される数字列を取り扱う問題。

解法
数式一発で解いても良かったけど、1001*1001と小さなサイズなので足し合わせて解いてみたり。
プロジェクトオイラーの問題は基本問題番号が大きいほど難しくなります。
このへんの低い問題はまだまだ中学生レベルです。
まあこのへんでも賢い解き方と私みたいなコンピュータパワーに頼った余り賢くない解き方に分かれるのでプログラムは結構奥深いものがあるようです。
備考
整数論の教育を受けたことはありません。
そのうえ整数論は苦手分野なため、プロジェクトオイラーの問題を解くときはマシンパワーに頼ったコードを書いてます。
貴方が私のサイトのコードを見てダサいコードだと感じたならそれは実際ダサいコードです。


 calc(501,Sum):-!,write(Sum).
 calc(N,Sum):-Sum1 is Sum+(2*N-1)*(2*N-1)*4+10*2*N,
 	N1 is N+1,
 	calc(N1,Sum1).
 main28:-calc(1,1).


*Problem 29 「個別のべき乗」 †
http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2029
a^bで2<=a<=100,2<=b<=100
a^bで重複するものを省くと何個の数ができるか。


解法
この問題で他の人の回答も見てみたらpythonがなかなかクールな言語だと感心。
気楽なリスト処理にstd::setみたいな便利機能が非常に使いやすそう、かっこいい言語だ。
コードを短く解くならhttp://sozorogusa.blogspot.jp/2012/07/problem29.htmlみたいなコードが一番賢い。
Prologで書いてもsortで消し去ってlengthで長さを測ればいいのでほとんど同じ処理になる。


 pow_list(N,N,[],100):-!.
 pow_list(N,ReN,[ReN|Result],C):-C1 is C+1,pow_list(N,Re,Result,C1),ReN is Re*N.
 calc(List):-
 	between(2,100,A),
 	pow_list(A,_,List,1).
 main29:-findall(List,calc(List),List1),
 	flatten(List1,Ans1),
 	sort(Ans1,Ans2),
 	length(Ans2,AnsLen),write(AnsLen).


もうひとつ、めんどくさい発想によるコードもあります。
この問題は2~100までの間にある素数を次元名とするベクトルの問題です。
a^n=b^m同じ数になるということはa,bを素因数分解して素数ベクトル空間にマップすれば、ベクトルa*n=ベクトルb*m、n<100、m<100となるということを指します。
つまり一バイト変数しかない言語でも配列さえ使えればそのまま解けるということを指します。
これを検討すれば問題は解けますがめんどくさいのでそんなコードは書きません。





*Problem 30 「各桁の5乗」 †
驚くべきことに, 各桁を4乗した数の和が元の数と一致する数は3つしかない.

1634 = 1^4 + 6^4 + 3^4 + 4^4
8208 = 8^4 + 2^4 + 0^4 + 8^4
9474 = 9^4 + 4^4 + 7^4 + 4^4
ただし, 1=1^4は含まないものとする. この数たちの和は 1634 + 8208 + 9474 = 19316 である.

各桁を5乗した数の和が元の数と一致するような数の総和を求めよ.

解法
コードを高速化して解いてみた。
体感速度的には十分な速度。
まあこんなものかな。
各桁の5乗和は密度が低いのでそちらの組み合わせだけを検証して速度アップ。


 check([_|_],0):-!,fail.
 check([],0):-!.
 check(Nums,Sum):-
 	T is Sum mod 10,
 	select(T,Nums,Rest),!,
 	Sum1 is Sum//10,
 	check(Rest,Sum1).
 search(_,_,_,Sum,_):-355000<Sum,!,fail.
 search(Nums,_,_,Sum,Sum):-
 	check(Nums,Sum).
 search(_,_,6,_,_):-!,fail.
 
 search(Nums,P,Deep,Sum,Result):-
 	between(P,9,N),
 	Sum1  is Sum+N*N*N*N*N,
 	Deep1 is Deep+1,
 	search([N|Nums],N,Deep1,Sum1,Result).
 
 sum([],Sum,Sum):-!.
 sum([X|Rest],Sum,Result):-
 	Sum1 is Sum+X,sum(Rest,Sum1,Result).
 
 main30_2:-
 	findall(Re,search([],0,0,0,Re),List),
 	sum(List,0,Ans),Ans1 is Ans-1,write(Ans1).