「prolog勉強プロジェクトオイラー61～70」の編集履歴（バックアップ）一覧に戻る
プロジェクトオイラーの問題を堀江伸一さんがPrologで解くページ
*問い61
この問題は探索でといても十分間に合いますが。
スタック操作だけで動的計画法で解くというチャレンジをしてみました。
どの部分の処理も線形時間か n log(n)で終わるようになっています。
appendを使ってるところだけが効率が悪いくらいです。
calcP3(N1):-between(1,200,N),N1 is (N*(N+1))//2,999<N1,N1<10000.
calcP4(N1):-between(1,100,N),N1 is N*N,999<N1,N1<10000.
calcP5(N1):-between(1,100,N),N1 is (N*(3*N-1))//2,999<N1,N1<10000.
calcP6(N1):-between(1,100,N),N1 is (N*(2*N-1)),999<N1,N1<10000.
calcP7(N1):-between(1,100,N),N1 is (N*(5*N-3))//2,999<N1,N1<10000.
calcP8(N1):-between(1,100,N),N1 is (N*(3*N-2)),999<N1,N1<10000.
mysweep([],[]):-!.
mysweep([[NNum,FNum,Sum]],[]):-integer(NNum),integer(FNum),integer(Sum),!.
mysweep([[NNum,FNum,_],[NNum,FNum,Sum]|Rest],Result):-!,
mysweep([[NNum,FNum,Sum]|Rest],Result).
mysweep([_,Set|Rest],[Set|Result]):-
!,
mysweep([Set|Rest],Result).
calc2([NN,FN,Sum],[[NN,NB]|_],[NB,FN,Sum1]):-
Sum1 is Sum+NN*100+NB.
calc2(Set,[_|Rest],Result):-
calc2(Set,Rest,Result).
calc([],_,Ans,Perm):-!,search(Ans,Perm).
calc(_,[],Ans,Perm):-!,search(Ans,Perm).
calc([[NN,FN,Sum]|Rest],[[NN,NB]|Rest1],Ans,Perm):-
!,
findall(Re1,calc2([NN,FN,Sum],[[NN,NB]|Rest1],Re1),List),
append(Ans,List,Ans1),
calc(Rest,[[NN,NB]|Rest1],Ans1,Perm).
calc([[NN,_,_]|Rest],[[NF,NB]|Rest1],Ans,Perm):-
NN<NF,!,
calc(Rest,[[NF,NB]|Rest1],Ans,Perm).
calc(List,[_|Rest],Ans,Perm):-
!,
calc(List,Rest,Ans,Perm).
check([[A,A,Ans]|_],Ans):-!.
check([_|Rest],Ans):-
check(Rest,Ans).
search([],[]):-!,fail.
search(Ans,[]):-check(Ans,Ans1),write([ansis,Ans1]).
search(Ans,Perm):-
sort(Ans,Ans1),
[Top|_]=Ans1,
mysweep(Ans1,Re1),
Ans2=[Top|Re1],
select(L,Perm,Rest1),
calc(Ans2,L,[],Rest1).
to_date3([],[]):-!.
to_date3([Y|Rest],[[N2,N1,Y]|Result]):-
N1 is Y // 100,
N2 is Y mod 100,
to_date3(Rest,Result).
to_date2([],[]):-!.
to_date2([Y|Rest],[[N2,N1]|Result]):-
N2 is Y // 100,
N1 is Y mod 100,
to_date2(Rest,Result).
to_date([],[]):-!.
to_date(List,Re2):-
member(L,List),
to_date2(L,Re1),
sort(Re1,Re2).
main61:-
findall(Y,calcP3(Y),L3),
findall(Y,calcP4(Y),L4),
findall(Y,calcP5(Y),L5),
findall(Y,calcP6(Y),L6),
findall(Y,calcP7(Y),L7),
findall(Y,calcP8(Y),L8),
findall(Re,to_date([L3,L4,L5,L6,L7],Re),Perm),
to_date3(L8,L88),
search(L88,Perm).
*問い62
Nを小さい数字から始めN^3の各桁の数字の出現数でカウントしていき大きな数字へ続けます。
条件を満たせばそれを暫定解として探索していき暫定解*10＜N＾３となればそれより小さい答えはありません。
ｃ++のstd::mapにあたる機能がピュアPrologにはないので線形探索をしておりコード実行速度は遅いです。
std::mapに当たる機能があったら非常に高速に解ける問題です。
mycount(0,Temp,Result):-sort(Temp,Result),!.
mycount(N,Temp,Result):-
N1 is N//10,
M is N mod 10,
(select([M,C],Temp,Rest) ->
C1 is C+1,
mycount(N1,[[M,C1]|Rest],Result);
mycount(N1,[[M,1]|Temp],Result)).
search(N,Ans,_):-
0<Ans,
T is Ans*10,
N3 is N*N*N,
T<N3,
!,
write(Ans).
search(N,Ans,Counts):-
N3 is N*N*N,
N1 is N+1,
mycount(N3,[],Re1),
select([Re1,C,Min],Counts,Rest),!,
(4 =< C ->
((Ans>Min;Ans=:=0) -> search(N1,Min,Counts);
search(N1,Ans,Counts));
C1 is C+1,
search(N1,Ans,[[Re1,C1,Min]|Rest])
).
search(N,Ans,Counts):-
N3 is N*N*N,
N1 is N+1,
mycount(N3,[],Re1),
search(N1,Ans,[[Re1,1,N3]|Counts]).
*問い63
この問題は10より大きな数を考えると10＾２＞2桁なのでそれより大きな数の乗数は考慮する必要はありません。
1～9まで検証すれば十分です。
Logを使った方がすっきりしますが桁数を簡単に数えてときました。
keta_size(0,0):-!.
keta_size(N,Result):-
N1 is N//10,
keta_size(N1,Re),
Result is Re+1.
search2(A,B,N):-N is A^B,keta_size(N,L),B>L,!,fail.
search2(A,B,N):-N is A^B,keta_size(N,L),B=:=L.
search2(A,B,N):-B1 is B+1,search2(A,B1,N).
search(N):-
between(1,9,A),
search2(A,1,N).
main63:-findall(N,search(N),List),sort(List,List1),length(List1,Len),write(Len).
*問い64
10000以下の平方数でない数を平方根を連分数展開で計算する。
これを連分数展開した時、連分数展開が奇数周期になる物の個数を答えよ。
gcd(0, B, B) :-!.
gcd(A, B, G) :- A \== 0, R is B mod A, gcd(R, A, G).
calc(N,NSqrt,U,DR,Deep,List,Result):-
D is N-DR*DR,
gcd(D,U,G),
D1 is D//G,
U1 is U//G,
L1 is floor((NSqrt+DR)*U1/D1),
L2 is L1*D1,
UR2 is L2-DR,
Deep1 is Deep+1,
(member([L1,UR2,Deep2,D1],List) ->
!,Result is (Deep1-Deep2) mod 2
;calc(N,NSqrt,D1,UR2,Deep1,[[L1,UR2,Deep1,D1]|List],Result)).
search(10001,Ans):-!,write(Ans).
search(N,Ans):-
sqrt(N,N1),
N=:= floor(N1)*floor(N1),
!,
N2 is N+1,
search(N2,Ans).
search(N,Ans):-
sqrt(N,NSqrt),
FloorSQ is floor(NSqrt),
calc(N,NSqrt,1,FloorSQ,0,[],Result),
Ans1 is Ans+Result,
N1 is N+1,
search(N1,Ans1).
*問い65
素朴に定義通り求めてみました。
gcd(0, B, B) :-!.
gcd(A, B, G) :- A \== 0, R is B mod A, gcd(R, A, G).
add(U,D,U1,D1,ReU,ReD):-
U2 is U*D1+U1*D,
D2 is D*D1,
gcd(U2,D2,G),
ReU is U2//G,
ReD is D2//G.
calc([X],U,D):-
integer(X),
U is 1,
D is X.
calc([X|Rest],U1,D1):-
calc(Rest,U,D),
add(U,D,X,1,D1,U1).
calc_list(M):-
between(0,98,N),
(1=:=N mod 3->
M is N //3*2+2;
M is 1).
keta_sum(0,0):-!.
keta_sum(N,Result):-
N1 is N//10,
keta_sum(N1,Re),
Result is Re+N mod 10.
main65:-
findall(M,calc_list(M),List),calc(List,U,D),add(U,D,2,1,U1,_),
keta_sum(U1,Ans),write(Ans).
*問い66 http://projecteuler.net/problem=66
この問題は色々式変形したり漸化式が立たないか考えたけど自力では解けなかった問題。
ネットに転がっていたjavascriptコードをＰｒｏｌｏｇコードに翻訳して100％カンニングして解いた。
驚くほど高速に答えが出てすこし驚いた。
この手法は連分数展開というものらしい。
大雑把な考え方までは理解したけれど、細かい計算は理解してない状態。
D=(ｘ＾２-1)/y^2はｙが大きくなると
x^2/y^2と見分けがつかなくなってくる。
√Ｄ＝x/yを満たす分数を求める問題で近似するという発想らしい。
聞けばなるほどの着想だが厳密解といいう考え方にこだわってると中々出てこない発想だな。
一つ勉強になった。
calc(_,_,_,_,E,F,K,N,_,P2,_,_,X,Result):-
B1 is (X-E*E)/F,
C1 is floor((K+E)/B1),
A1 is C1,
A1=:=2*K,0=:=N mod 2,!,Result is P2.
calc(_,_,_,_,E,F,K,N,P1,P2,Q1,Q2,X,Result):-
B1 is (X-E*E)/F,
C1 is floor((K+E)/B1),
D1 is C1*B1-E,
A1 is C1,
E1 is D1,
F1 is B1,
P3 is A1*P2+P1,
Q3 is A1*Q2+Q1,
N1 is N+1,
calc(A1,B1,C1,D1,E1,F1,K,N1,P2,P3,Q2,Q3,X,Result).
calc_first([Result,D]):-
between(2,1000,D),
sqrt(D,T),
K is floor(T),
K2 is K*K,
D\==K2,
A is K,
E is K,
F is 1,
N is 1,
P1 is 1,
P2 is K,
Q1 is 0,
Q2 is 1,
calc(A,_,_,_,E,F,K,N,P1,P2,Q1,Q2,D,Result).
main66_2:-
findall([D,X],calc_first([D,X]),List),
sort(List,List1),write(List1).
*問い67
問題自体は簡単。
ただPrologには破壊的代入がないせいで少しめんどくさい。
とりあえずProlog風味なファイル読み込みとは何かと考えてコーディングをしてみました。
データにゴミや入力ミスが入ってないから練習問題って楽でいいですね。
one_read(Result,State):-get0(C),((C =:= -1 ; C=:=10;C=:=32)->Result=[],State is C;
one_read(Re,State),Result=[C|Re]).
row_read([Re1|Result],ReState):-
one_read(Re,State),(Re==[]->Re1=[];to_num(Re,0,Re1)),((State=:= -1;State=:=10)->ReState is State,Result=[];
row_read(Result,ReState)).
all_read([Result|Re]):-row_read(Result,State),(State=:= -1 ->Re=[];
all_read(Re)).
to_num([],Sum,Sum):-!.
to_num([X|Rest],Sum,Result):-
Sum1 is Sum*10+(X-48),to_num(Rest,Sum1,Result).
calc([],Sums):-!,write(Sums).
calc([_|Rest],Sums):-
calc(Rest,Sums).
main67:-seen,see('e67.txt'),all_read(Re),write(Re),calc(Re,0).
