小屋畑ソフトウェアサービス ブログページ

ブログ

  • Project Euler Problem 52

    Problem52

    ECLiPSeで解きました。
    処理時間は自分のマシン(Intel(R) Core(TM) i5-3340M CPU @ 270GHz)で6.5秒くらい。
    見どころ？はECLiPSeのlogical-loopで書いた順列permutationと、propiaで数値を入れ替えた制約を作っているところ
    logical loop はループ用のサブ述語みたいなアホっぽい定義がなくなってコードがエレガントになるのでSWI-Prologでもぜひ正式採用してほしい。
    桁のところは一般化してなく、1桁、2桁、3桁…　と手入力で増やしながら確認しました
    （bb_minでそれぞれの桁での最小の数値が取れるようにしてます）。

    プログラム：

    :-lib(ic).
:-lib(propia).
:-lib(branch_and_bound).

go(Digits,[X1num,X2num,X3num,X4num,X5num,X6num]) :-
	length(X1Lst,Digits),
	X1Lst#::[0..9],
	lst2num(X1Lst,X1num),
	X1num #> 10 ^ (Digits-1),
	X2num #= X1num * 2,
	X3num #= X1num * 3,
	X4num #= X1num * 4,
	X5num #= X1num * 5,
	X6num #= X1num * 6,
	permutation(X1Lst,X2Lst) infers most,
	lst2num(X2Lst,X2num),
	permutation(X1Lst,X3Lst) infers most,
	lst2num(X3Lst,X3num),
	permutation(X1Lst,X4Lst) infers most,
	lst2num(X4Lst,X4num),
	permutation(X1Lst,X5Lst) infers most,
	lst2num(X5Lst,X5num),
	permutation(X1Lst,X6Lst) infers most,
	lst2num(X6Lst,X6num),
	labeling([X1num,X2num,X3num,X4num,X5num,X6num]).
	
lst2num(L,N):-
	(foreach(Num,L),fromto(0,In,Out,N),count(I,0,_)
		do
		Out #= In + 10^I * Num
	).

permutation(A,B):-
	(fromto(A,AIn,AOut,[]),fromto(B,[BFirst|BRest],BRest,[])
		do
		select(BFirst,AIn,AOut)
	).

    実行結果（解答伏せます）：

    ?- bb_min(go(XXXXXXX, [N1, N2, N3, N4, N5, N6]), N1, Option).
N1 = XXXXXXXX
N2 = XXXXXXXX
N3 = XXXXXXXX
N4 = XXXXXXXX
N5 = XXXXXXXX
N6 = XXXXXXXX
Option = bb_options(continue, -1.0Inf, 1.0Inf, 1, 1, 0, 0, one, _466, _467, _468)
Yes (6.59s cpu)

  • Project Euler Problem79

    久しぶりにProject Eulerやってみました。
    制約論理プログラミングで解きました。ほとんど頭使ってない…
    SWI-Prologです。

    Problem79

    :-use_module(library(clpfd)).

solve(A,Len):-
	length(A,Len),
	A ins 0..9,
	choice(A,[3,1,9],IdxLst1,Len),
	choice(A,[6,8,0],IdxLst2,Len),
	choice(A,[1,8,0],IdxLst3,Len),
	choice(A,[6,9,0],IdxLst4,Len),
	choice(A,[1,2,9],IdxLst5,Len),
	choice(A,[6,2,0],IdxLst6,Len),
	choice(A,[7,6,2],IdxLst7,Len),
	choice(A,[6,8,9],IdxLst8,Len),
	choice(A,[7,6,2],IdxLst9,Len),
	choice(A,[3,1,8],IdxLst10,Len),
	choice(A,[3,6,8],IdxLst11,Len),
	choice(A,[7,1,0],IdxLst12,Len),
	choice(A,[7,2,0],IdxLst13,Len),
	choice(A,[7,1,0],IdxLst14,Len),
	choice(A,[6,2,9],IdxLst15,Len),
	choice(A,[1,6,8],IdxLst16,Len),
	choice(A,[1,6,0],IdxLst17,Len),
	choice(A,[6,8,9],IdxLst18,Len),
	choice(A,[7,1,6],IdxLst19,Len),
	choice(A,[7,3,1],IdxLst20,Len),
	choice(A,[7,3,6],IdxLst21,Len),
	choice(A,[7,2,9],IdxLst22,Len),
	choice(A,[3,1,6],IdxLst23,Len),
	choice(A,[7,2,9],IdxLst24,Len),
	choice(A,[7,2,9],IdxLst25,Len),
	choice(A,[7,1,0],IdxLst26,Len),
	choice(A,[7,6,9],IdxLst27,Len),
	choice(A,[2,9,0],IdxLst28,Len),
	choice(A,[7,1,9],IdxLst29,Len),
	choice(A,[6,8,0],IdxLst30,Len),
	choice(A,[3,1,8],IdxLst31,Len),
	choice(A,[3,8,9],IdxLst32,Len),
	choice(A,[1,6,2],IdxLst33,Len),
	choice(A,[2,8,9],IdxLst34,Len),
	choice(A,[1,6,2],IdxLst35,Len),
	choice(A,[7,1,8],IdxLst36,Len),
	choice(A,[7,2,9],IdxLst37,Len),
	choice(A,[3,1,9],IdxLst38,Len),
	choice(A,[7,9,0],IdxLst39,Len),
	choice(A,[6,8,0],IdxLst40,Len),
	choice(A,[8,9,0],IdxLst41,Len),
	choice(A,[3,6,2],IdxLst42,Len),
	choice(A,[3,1,9],IdxLst43,Len),
	choice(A,[7,6,0],IdxLst44,Len),
	choice(A,[3,1,6],IdxLst45,Len),
	choice(A,[7,2,9],IdxLst46,Len),
	choice(A,[3,8,0],IdxLst47,Len),
	choice(A,[3,1,9],IdxLst48,Len),
	choice(A,[7,2,8],IdxLst49,Len),
	choice(A,[7,1,6],IdxLst50,Len),
	label(A).


choice(A,[P1,P2,P3],[N1,N2,N3],Cnt):-
	[N1,N2,N3] ins 1..Cnt,
	element(N1,A,P1),
	element(N2,A,P2),
	element(N3,A,P3),
	N1 #< N2,
	N2 #< N3.

    回答伏せます

  • ドローン

    なんか流れでドローンの制御をするソフトに携わることになってドローンの事情など何も知らない状態から苦労しながらどうにかこうにか設計にチャレンジしています。

    今のドローンはGPSと画像処理などを駆使してムチャクチャ頭良い処理ができるようになっていて、例えばコントローラーのタブレットで走る車などにターゲット指定したら移動するターゲットを常時距離を保ったまま追いかけ続け、カメラ（ジンバル）もターゲットに常時固定しながら動画を撮れる。

    素早く動くバイクなどが常時カメラに納まり続けているの動画を見てると、例えばもしこれが人を追いかけている状況だとしたらと考えると恐ろしい技術で、ジンバルに拳銃を固定したとして、遠隔で引き金を引くなんてプログラムは簡単にできるのだから、簡単にターミネーターのような殺人機械ができそうです。

    POIとかいうシステムで仏像とかのターゲットを指定したら、その周りをドローンがくるくる撮影しながらながら回り続けて、リアルな３ｄのモデリングデータ作成も出来る。

    ドローンは姿勢や位置の制御のために各種のセンサが搭載されているのだが、どうもこれらが安価に手に入るようになったことが安くドローンを作ることができるようになった理由らしくて、これらのセンサが安価につくれるようになった理由はスマホが大量生産されるようになったかららしい。

    確かにドローンに搭載されているGPSセンサ、気圧センサ、加速度センサ、磁気センサなどはどれもスマホに搭載されている。

    シヴィライゼーションというゲームで技術の順番を記述したテクノロジーツリーという考え方があり、例えば「法律」という技術は「アルファベット（文字）」という技術を取得していないと開発できないように順番があるのですが、さしずめ「ドローン」という技術の前提技術として「スマートフォン」が来る感じであろう。ちなみに「アルファベット」の前提技術は「筆記」となっている。

    20年前には自分がドローン制御のプログラムを書くことになるとは思いもしなかった。
    時代の進み方の速さにめまいがしそうだし、自分も末端の人間としてそのような進歩にほんの少しだけ関わっている

    勉強していてもオープンソース系でそのような技術のソフトウェア部分が構成されているのが多い印象で、技術進歩のスピードが無茶苦茶はやいのにはオープンソースコミュニティが1役も2役も買っていると思う

  • [Prolog] DCGで双方向性のある述語を作成

    最近stackoverflowでPrologの質問に良く回答しているのですが、面白そうな問題があり、自分でもプログラムを作ってみました。

    質問

    内容は、[1,2,3,4,1,4,7,5,7,8,9] のようなリストが与えられた際に「続いている数字は[最初の数字,最後の数字]の形で表し、そうでないものは単体の数字で表したようなリストを生成する　という問題です。
    [1,2,3,4,1,4,7,5,7,8,9]は[[1,4],1,4,7,5,[7,9]] に変換されます。

    作っているうちにどんどんアイデアがわいてきて、はじめは長いプログラムだったのがDCG使ってエレガントな方法で実装する方法を思いつき、最後には入出力逆転してもOKなような非常に汎用的な形で作成できました。

    ECLiPSeで作成(多分他の環境でも問題なく動きます）:

    chunk_list([])-->[],!.
chunk_list([Chunk|Rest])-->chunk(Chunk),!,chunk_list(Rest).

% sequence
chunk([First,Last])-->sequence([ First , Last ] ),{ First < Last }. 

% isolated number
chunk(Val)-->sequence([Val,Val]).               

sequence([First,Last]) --> 
                 { ( number ( First ) , number( Last ) ) -> First < Last ; true } ,
                 [ First ] ,
                 { succ ( First , Second ) } , sequence( [ Second , Last ] ) .

sequence([Val,Val])-->[Val].

    実行結果：

    [eclipse 4]: phrase(chunk_list([2,[2,6],[3,5],3,7,[1,3]]),Ret).

Ret = [2, 2, 3, 4, 5, 6, 3, 4, 5, 3, 7, 1, 2, 3]
Yes (0.00s cpu)

[eclipse 5]: phrase(chunk_list(Ret),[3,4,5,2,1,6,7,88,9,4,5,6,7,2,1,2,3]).

Ret = [[3, 5], 2, 1, [6, 7], 88, 9, [4, 7], 2, [1, 3]]
Yes (0.00s cpu)

[eclipse 6]: phrase(chunk_list([[2,4],A,[2,8],3]),[2,B,4,6,2,C,4,5,6,7,D,E]).

A = 6
B = 3
C = 3
D = 8
E = 3
Yes (0.00s cpu)

    最後のマッチング例は「やってみたら、できた」という実験で、自分の予想を上回る汎用性が実現できていてびっくりしました。

    入出力がどちらもリストの述語を作成するときは、工夫すると双方向性が実現できることがあると思います。

    Prologの真髄に触れたような気がしました。

  • [Prolog]最近の勉強

    最近あんまりPrologのブログ書いてませんが、ずっと研究してます(具体的にブログに書けるような結果が出てない）現在はECLiPse（IDEのほうじゃなくて制約プログラミングのほう）の勉強してます。

    書くネタないので今までのPrologの勉強の流れを説明すると以下のような感じで進んできました(スパンは数年単位）

    広井さんのページでPrologの基礎を把握しいくつかのパズルのソルバーを作る（20年近く前に大学の授業で受けていたがあまり深いところまで理解していなかった）

    ↓

    SWI-Prologのマニュアル読み込むうちCLPFDを発見し衝撃を受け、勉強開始

    ↓

    ProjectEulerでPrologで問題を解くうち、正答者スレでECLiPseを用いてSWIの自分のプログラムよりはるかに速いスピード(数十倍）で結果を出している人を発見して衝撃を受ける。プログラムの内容自体はCLPFDを使用した自分のとそんなに変わらなかった。

    ↓

    ECLiPseの勉強を始める。だんだんSWI-PrologのCLPFDがSicstusやECLiPseの制約ソルバーに比べると「オモチャ」レベルだったということがわかってくる。また、SWI-PrologのCLPFD開発者のmarkus triskaさんは開発の継続をsicstusで行うみたいな記述を最新のSWIのマニュアルで見つけた。

    また、ECLiPse は node.jsを使用してjavascriptで呼び出せるらしく、ECLiPseの勉強はWeb系の仕事でのちのち面白い威力を発揮するかもしれないと期待を抱く。

    ↓

    ECLiPseには、自分がSWIのCLPFDでどうしてもできなかった自作制約の作成がPropiaとCHRの2種類の方法で可能なことを知り、そちらの勉強を開始。SWIのCLPFDの自作制約の作成方法はマニュアルが意味不明だった（まだ完成してないとも書いてある）

    ↓

    Propiaのほうは大体の概要を把握。任意の述語をそのまま制約にできて、Propagationの厳密さをinferという述語の引数でコントロールするという仕組みらしい

    ↓

    CHRの勉強に移る。ルール自体は少ないが結構難しくてどのように作ればよいかなかなか理解できない。CHRはHaskellやJavaなどでも使用できるようだ。

    ↓

    CHRが並列計算に向いていると知り、ハードウェアで実行したりGPUを使用したりして解くといった論文をネットで見つけて読み始め、少しECLiPseの勉強から浮気する。

    ↓

    「CHRをGPUを用いて計算する｣方法を記載した論文を見つけるが実際にやった時の速度比較が載ったものを見つけることができなかったので、そういう論文が出ないかを定期的にウォッチすることに決め、とりあえず論文リサーチはやめる。CHRをハードウェアでやった論文では数万倍とかのスピードで解いてるの見つけたが良く読むとCHRの機構だけではなくプログラムのほうもハードウェア化していてちょっとちがうだろという感じだった（論文で使用している最大公約数求めるみたいな単純なロジックしかハードウェア化できなくて現実的ではないし、そもそもそれってCHR自体の速さではない）

    ↓

    またECLiPseのCHRの勉強に戻る

    △いまここです。

    まあ日本語の情報がないです。
    ずーーっと英語読んでます。ECLiPse、日本で使ってる人100人もいないんじゃなかろうか。日本のECLiPse界では望まずして第一人者になりそうです。もしそういう人たちいれば仲間に入れてほしい。誰か一緒に勉強しませんか？情報交換して一緒に研究の成果を教えあいましょう

    CHRはルール自体は3種類しかなくて単純なんだけど以下が良くわかってなくてまだ自作の制約を作るまで至ってない

    ・labelingの振る舞い
    ・制約の述語の変数に「値を設定した時？」と「新しい制約がセットされたとき？」に再度bodyがfireされるらしいがそこらへんの振る舞い

    多分CHRの制約の形が｢自分がCLPFDを勉強してきた上で持っている制約のイメージ」とかなりずれているのが理解が滞ってる理由だと思います。

    自分の数年の目標であるTantrixというパズルのソルバーをECLiPseの自作制約でいつか作れれば良いなあ。手続き型だと現在のスキルで多分作れるが制約プログラミングで作るという目標がありわざとそういう手法では作ってない

    …………………………うーん。

    この記事書いてるうちになんで自分がCHRを必死に勉強してるのか分かんなくなってきた…
    Propiaで間に合ってるような気がしてきたので、そっちでやるかな…。CHRは「つかうとすごいスピードが出るケースが見つかる(並列処理とかで）」とか「CHRで書かれためちゃくちゃエレガントなコードを見つける」みたいな状況になってから改めて勉強しようかな。
    CHRの良さをご存じの方いらっしゃればぜひ koyahata*koyahatataku.com(*を@に変更) にメール下さい。そうでなくてもProlog関連で勉強中の人はぜひ連絡取り合って一緒に勉強しましょう

    ブログまとめるのも自分の行動を見直す良い機会になりますね

    よし方向転換

  • C#のWPD(MTP)操作

    C#使ってMTPデバイスに読み書きするプログラム作ってるんですがこれがまあめんどい。

    WPDというAPIを使うんですがWeb上に全然動くサンプルがないし、Microsoftのページもかなり不親切な構成になっていて
    簡単なことやろうとしても非常に苦労します。

    極めつけはGetDevicesという接続デバイスを列挙する関数にバグがあって１デバイスしか取得できず、それを解決するためには
    Dllを逆アセンブル→修正→再アセンブルが必要という始末…
    参照

    フォルダ移動・アップロード・ダウンロード・削除といった基本機能ができるようになったプロジェクトをアップしますので、同じことやってる人はこれ使って時間無駄にしないようにしてください…いやーしょうもないことで時間使ってしまった。
    サンプルプロジェクト
    この程度のソースが全然手に入んないんです…

    Visual Studio 2012プロジェクトです。

  • PrologでForループ

    SicstusやECLiPSeではイテレータとかいうのを使ってDo-Loopが記述できるのですが、もともとECLiPSeで実現されていた論理プログラミングのループ(Logical Loop)を紹介する以下の論文を2002年にSchimpfさんという人が書いていて、これでLogical Loopが広まったようです？

    Logical Loop

    まだ深く読み込んでませんが、SWIにもなんかの方法で実装できるような気もする… lambdaライブラリとかを使う？

  • 発見

    CLP(B)を使ってProjectEuler の No.172を解いている記事を見つけた。

    CLPFD の開発者のMarkus Triskaさんが書いている。

    The Boolean Constraint Solver of SWI-Prolog: System Description

    他にも No.161 に関する記載もありますね。

    …しかし律儀に全数を数え上げる方法を使ってるからムチャクチャ効率悪くて、メモリ20Gbyte積んだマシン使って7時間もかかって解いてる…ダメだなこりゃ…

  • Probrem 185 Number Mind

    Project Eulerの問題はものすごくたくさんあるし（2017.05.12現在で600問弱）自分の目的はProlog関連技術の威力を試すためにやっている側面が強く全問制覇が目標ではないので、これからはまずそれに適した問題から解いてゆこうと思う。

    Probrem 185 Number Mind

    この問題は「ザ・制約論理プログラミング用問題」というくらい制約論理プログラミングに向いていると感じ手を付けた。ほとんど探索の工夫をせずに問題の定義のみ書いて、ちゃんと解けている。

    今回プログラムを作っていてハマってしまった箇所があり、勉強になった。
    ・術語内で制約論理用変数を使いそれを呼び元に返さない状態でlabelingを行うと、それらの変数がガベージコレクトされておかしくなる。
    ・上記のエラーメッセージはguess述語の呼び出しをかなり減らした状態から順々に増やしてゆかないと見ることができない（全て書いた状態で実行すると一晩返ってこないプログラムとなり時間をかけて正しく動いているようにも見え、プログラムが間違っているという手がかりは全く得られなかった）

    上記のような手続き型言語と異なる独自の開発ノウハウが結構あるので、これからも要勉強である。

    labelingのオプションでも相当速度が変わる。基本はffもしくはffcをつけるで良いと思う。

    この問題は Difficulty Rating が55% で、 2017年5月12日現在で解いた人が 2273 人しかいない問題だそうで、今まで自分が解いた中では一番難しい問題のようだ。

    個人的には今までチャレンジかつ正答した問題の中では 191 Prize String が一番難しく、これは1週間位組み合わせ・順列系のサイトを調べまくって解いた。自分の解き方はエレガントでなかったが、正答者のみ見れるスレッドによると動的計画法でエレガントかつ超スピードで解けるようだ。

    Number Mindの問題は前述のハマり以外はとくに問題はなかった。正答者スレによるとECLiPSeとかいう（IDEのほうではない）同じく制約論理ライブラリを使って解いている人がいた。その人のコードもほぼ問題の定義のみが記述されているエレガントなものだった。自分のスピードより数十倍速いようだったので、Eclipse恐るべしと思った。いつか勉強するかも。

    正答者スレを見ることにより、自分より遥かに優秀な取り組み方をしている人の使用言語やソースコードを確認できるのも Project Euler の良いところだと思う。

    ソースコード：

    :-use_module(library(clpfd)).
solve(Num):-
	length(Num,16),
	Num ins 0..9,
	
	guess(Num,2321386104303845,0,Hits00),
	guess(Num,6375711915077050,1,Hits01),
	guess(Num,6913859173121360,1,Hits02),
	guess(Num,3847439647293047,1,Hits03),
	guess(Num,3174248439465858,1,Hits04),
	guess(Num,8157356344118483,1,Hits05),
	guess(Num,4895722652190306,1,Hits06),
	guess(Num,4513559094146117,2,Hits07),
	guess(Num,5616185650518293,2,Hits08),
	guess(Num,2615250744386899,2,Hits09),
	guess(Num,6442889055042768,2,Hits10),
	guess(Num,2326509471271448,2,Hits11),
	guess(Num,5251583379644322,2,Hits12),
	guess(Num,2659862637316867,2,Hits13),
	guess(Num,5855462940810587,3,Hits14),
	guess(Num,9742855507068353,3,Hits15),
	guess(Num,4296849643607543,3,Hits16),
	guess(Num,7890971548908067,3,Hits17),
	guess(Num,8690095851526254,3,Hits18),
	guess(Num,1748270476758276,3,Hits19),
	guess(Num,3041631117224635,3,Hits20),
	guess(Num,1841236454324589,3,Hits21),
	
	Hits=[	Hits00,Hits01,Hits02,Hits03,Hits04,
			Hits05,Hits06,Hits07,Hits08,Hits09,
			Hits10,Hits11,Hits12,Hits13,Hits14,
			Hits15,Hits16,Hits17,Hits18,Hits19,
			Hits20,Hits21],
	
	flatten([Num,Hits],Flatten),
	labeling([ff,bisect],Flatten).
	
guess(Num,CompareNum,HitNum,Hits):-
	number_codes(CompareNum,Codes),
	maplist(code_num,Codes,Nums),
	%write(Nums),nl,
	maplist(chk_hit,Num,Nums,Hits),
	sum(Hits,#=,HitNum).

chk_hit(A,B,H):-
	A#=B #<==> H.
	
code_num(Cd,Num):-
	Num is Cd - 48.

    実行結果：

    ?- time(solve(Num)),write(Num).
% 585,922,243 inferences, 43.774 CPU in 43.975 seconds (100% CPU, 13385202 Lips)
**********解答伏せてます**********
% 176,962,326 inferences, 13.057 CPU in 13.140 seconds (99% CPU, 13552767 Lips)
false.

  • Problem2　偶数のフィボナッチ数

    問題

    % problem 2
solve:-
	fibo([2,1],Fibo),
	even_sum(Fibo,Sum),
	write(Sum).
	
even_sum([],0):-!.
even_sum([First|Rest],Sum):-
	First mod 2 =:= 0,
	!,
	even_sum(Rest,Sum1),
	Sum is First + Sum1.
even_sum([_|Rest],Sum):-
	even_sum(Rest,Sum).
	
	
fibo([First|Rest],Rest):-
	First >= 4000000,
	!.
fibo([P1,P2|Rest],Out):-
	P0 is P1 + P2,
	fibo([P0,P1,P2|Rest],Out).

    実行結果：
    1 ?- solve.
    ***********************(解答伏せます)**********************
    true.