■油売り算

そんなことして、油売ってる場合ではない

■Deriving Backtracking Monad Transformers

http://citeseer.ist.psu.edu/303981.html

を読んだ。Backtracking~とかはどうでもよくて、論文でTerm implementationとContext passin implementationというのが出てくるけど、この二つの定義の関係がよくわからん。

Term implementationって、例えば、GADTが使えれば、State Monadは次のようにも定義できるよね、とか、そういう話だと思う

data State' s a where
 Return :: a -> State' s a
 Bind :: forall b.(State' s b) -> (b->State' s a) -> (State' s a)
 Get :: State' s s
 Set :: s->State' s ()

GADT使わないTerm implementationは、例えば

{-# OPTIONS -fglasgow-exts #-}
--code from http://www.haskell.org/pipermail/haskell-cafe/2005-January/008258.html
import Control.Monad
import Control.Monad.State hiding (runState)

data State' s a = Return a
                | Get (s->State' s a)
                | Set s (State' s a)

instance Monad (State' s) where
  return = Return
  (Return x)>>=k = k x
  (Set s m)>>=k = Set s (m>>=k)
  (Get f)>>=k = Get (\s-> f s>>=k)

instance MonadState s (State' s) where
 get = Get Return
 put s = Set s (Return ())

runState :: State' s a -> s -> (a,s)
runState (Return a) s = (a,s)
runState (Set s1 k) s = runState k s1
runState (Get f) s = runState (f s) s

まあ、それはいいけど、上のように定義した型(っつーかモナド)と普通のStateモナドの定義の関係がよくわからん。

■Prolog User Golfing System

HaskellでPrologっぽいものを作ってみた。正しいPrologの挙動が分からないのでアレだけども

{-# OPTIONS -fglasgow-exts #-}
{-# OPTIONS -fno-monomorphism-restriction #-}
module Prolog where
import List(nub)
import Monad(liftM)
import Control.Monad.State
data Term = LVar String
          | Atom String
          | IntT Integer
          | Nil
          | Cons Term Term deriving (Eq,Show)

data Predicate = Eq Term Term
               | And Predicate Predicate
               | Or Predicate Predicate deriving Show

type Subst = [(String , Term)]

--Substの適用
subst :: Subst -> Term -> Term
subst [] t = t
subst (a:as) t = subst as (subst' a t)
 where
  subst' (n , t) (LVar v) = if (n==v) then t else (LVar v)
  subst' s (Cons x y)  = Cons (subst' s x) (subst' s y)
  subst' s x = x

--Substの合成
substC :: Subst -> Subst -> Maybe Subst
substC theta1 theta2 = substR (theta1++theta2)
 where
  substR theta = substR' [] (Just theta)
  substR' _ Nothing = 
   fail "failed to compose"
  substR' prev (Just theta) = 
   if (prev==theta) then (return theta) else (substR' theta $ next theta)
  next theta = 
   let n=[(n,subst theta t)|(n,t)<-theta] in 
    if (contradict n) then Nothing else (Just n)
  contradict theta = 
   all (\x->length x/=1) (nub[map snd $ filter (\(v,t)->v==v') theta|(v',t')<-theta])


occurs :: String -> Term -> Bool
occurs s (LVar n) = (s==n)
occurs s (Cons t1 t2) = (occurs s t1) || (occurs s t2)
occurs _ _ = False

unify ::  Term -> Term -> Maybe Subst
unify t1 t2 = unifyT [] t1 t2

unifyT :: Subst -> Term -> Term -> Maybe Subst
unifyT theta (LVar v1) (LVar v2) = 
 if (v1==v2) then (return theta) else (fail "can't unify")
unifyT theta (LVar v) t = 
 if (occurs v t) then (fail "can't unify") else (substC [(v,t)] theta)
unifyT theta t (LVar v) = 
 if (occurs v t) then (fail "can't unify") else (substC [(v,t)] theta)
unifyT theta (Cons t1 t2) (Cons t1' t2') =
 do eta<-unifyT theta (subst theta t1) (subst theta t1')
    unifyT eta t2 t2'
unifyT theta x y = 
 if (x==y) then (return theta) else (fail "can't unify")


solve' ::  Subst -> Predicate -> Maybe [Subst]
solve' init (Eq p q) = 
 do a<-unify (subst init p) (subst init q)
    ret<-substC init a
    return [ret]
solve' init (Or p q) =
 case (solve' init p) of
  Just ps -> case (solve' init q) of {Just qs->return (ps++qs);Nothing->return ps}
  Nothing -> solve' init q
solve' init (And p q) = 
 case (solve' init p) of
  Nothing -> Nothing
  Just ps -> let tmp = (map (\x->solve' x q) ps) in
   if (all (Nothing==) tmp) then Nothing else return $ concat[x|Just x<-tmp]


solve :: (Monad m)=> Predicate -> m [Subst]
solve p = case (solve' [] p) of{Just x->return x;Nothing->fail "failed"}

Termの方は、整数とか、シンボルとか、リストを扱える。で、Predicate。例えば、

like(john , taro).
like(masaki , abeTakakazu).

は

like :: Term -> Term -> Predicate
like a b = (Or (And (Eq a (Atom "John")) (Eq b (Atom "taro")))
                 (And (Eq a (Atom "masaki")) (Eq b (Atom "abeTakakazu"))))

とかいう感じの式になる。で、

*Main> solve $ like (LVar "x") (LVar "y")
[[("x",Atom "John"),("y",Atom "taro")],
 [("x",Atom "masaki"),("y",Atom "abeTakakazu")]]

とか、そんな感じ。

リストのappendとかは、

append :: Term -> Term -> Term -> Predicate
append ps qs rs = (Or (And (Eq ps Nil) 
                           (Eq qs rs))
                      (And (And (Eq ps (Cons (LVar "X") (LVar "Xs"))) 
                                (Eq rs (Cons (LVar "X") (LVar "Ys"))))
                            (append (LVar "Xs") qs (LVar "Ys"))))

とかだと多分まずい。っつーのは、appendを実際に評価したとき、何もしないと、appendの中にあるappendが再帰的にどんどん展開されて、無限に長い論理式ができるけど、そんときに、(LVar "X")とか(LVar "Ys")が違うスコープで何回もでてきて、衝突するから。なんで、gensymみたいな関数で、毎回、違う名前を生成しないといけない。とりあえず、お手軽に解決。

gensym :: State Int Term
gensym = do{n<-get;put $ n+1;return $ LVar(show n)}

app :: Term -> Term -> Term -> State Int Predicate
app ps qs rs = 
 do{x<-gensym;
    xs<-gensym;
    ys<-gensym;
    n<-app xs qs ys;
    return  (Or (And (Eq ps Nil) (Eq qs rs)) 
                (And (And (Eq ps (Cons x xs)) (Eq rs (Cons x ys))) n))}

で、

とか。一応正しく答えがでてるような気がする。あと、

solve (Eq (IntT 4) (IntT 3))

はfailして、

solve (Eq (IntT 3) (IntT 3))

は、[]が返ってくる。

いくつか、問題があって、一つは異なる論理変数XとYの単一化を無条件で失敗にしている点。これは、色々ダメな気がするけど、どうすれば綺麗に解決できるのか分からなかった。もう一点は、(And p q)解くときに、まずpを解いて、次にqを解いてるけど、qがfailする場合にはfailとしちゃっていいんだけど、pの計算が終了しない場合、全体として計算が停止しなくなる。要するに、pとqの計算は、本来並列に実行されるべきなのに、逐次計算になってるのがまずい。と思ったんだけど、よく考えると、Prologも前から順に計算していくなぁ

SWI-Prologで挙動確認してみると

hoge(X,Y) :- X=Y,Y=X,X=3.

は正しい答えを求めてくれる。appendの方は

append(X,Y,Z) :- X=[],Y=Z.
append(X,Y,Z) :- append(X1,Y,Z1),X=[W|X1],Z=[W|Z1].

と

append(X,Y,Z) :- X=[],Y=Z.
append(X,Y,Z) :- X=[W|X1],append(X1,Y,Z1),Z=[W|Z1].

は理屈上は同じ振る舞いをすべきだという気がするけど、前者で

?- append([],X,[1,2,3,4])

とかやると、解を一個求めて、更に計算を続行しようとすると、処理が戻ってこない。まあ、実行順序不定にしてしまうと

?- print(abe),print(takakazu)

とかも、出力がおかしくなるわけだけど、なんだかなぁ

人は、何故いとも簡単に生活ペースが乱れてしまうのか

■量子回路

すーぱーおぺれーたー続き。前回適当に書きなぐったことは、とりあえず9割くらいは間違ってなかった。

superoperatorについては、

http://www.lorentz.leidenuniv.nl/quantumcomputers/literature/preskill_1_to_6.pdf

あたりを読むと詳しく書いてあるけども、

・なんで単なるpositive mapでなく、completely positive mapなのか

物理的には、系Aと独立な系Bがあって、系Bの密度行列が時間発展しないとき、自然に定まる系A+Bの時間発展は、やはりpositive mapに決まってる。完全正値写像になる理由はそんだけ。数学的に見ると、positive map同士のテンソル積は、もはやpositive mapになるとは限らないので、都合悪いとか、operator sum representationという形式でsuperoperatorを表すためには、completely positive mapであることが必要とか、そんなのが単なるpositive mapだとダメな理由。

・なんでトレースが保存するという条件ではなく、トレースが増加しないという、より緩い条件なのか?

上にあげた文献では、トレースが保存するという条件になってるのだった。Selingerという人のTowards a Quantum Programming Languageという論文読むと、次のようなことが書いてある

しかし、トレースが保存するという条件を緩めることと、停止しない計算を扱えることに、どういう繋がりがあるのか分からない。

まあ、とりあえず、そういう感じで8割くらいは納得できた。

Haskellでの実装 in Sabryの論文。下準備編

type Vec a = a -> Complex Double

は、型aに属する値の線型結合で張られる複素ベクトル空間を表す。イメージは

Vec Bool = {a|True> + b|False> | a , b is Complex Numer}

とか、そんな感じ。Vecは殆どモナドだけども、一般には(>>=)が定義できない。これを定義するためには、Vec aの基底の集合を得る必要があって、それは、型aに属する値の全体なんだけども、Haskellではそれを取得する術はない。仕方ないので、

class Eq a => Basis a where

basis :: [a]

とか基底の全体を取得できる型クラスを定義する。BoolはBasisのinstanceで、単に

instance Basis Bool where

basis = [True , False]

で、「伝統的には」

type Lin a b = a -> Vec b

を使ってたけど、

type Super a b = (a,a) -> Vec (b,b)

使った方がよくね?という主旨らしい。Lin a bはVec aからVec bの線型写像全体で、Vec a->Vec bにしない理由は、そうすると線形でない関数まで入ってきてしまうから。線型写像は基底の行き先を決めれば全部の行き先が決まるので、a->Vec bでよい。同様に、Superの方は、意図はVec(a,a)からVec(b,b)への線型写像の全体。まあ、density matrixの棲息する空間は、波動関数の棲息する空間を二つテンソルしたものじゃなくて、波動関数の棲息する空間上のendomorphismの全体、要するに行列環と考えるべきだと思うけど、そのへんは、まあいい

さて、SuperはArrowになるよ！って書いてあるけど、Vecが厳密にはモナドでないのと同じ理由で、厳密にはArrowではないし、SuperがArrowであるという意味では、LinもArrowになるんだけど、そこは突っ込んではいけないの世界。論文では、fun2linとか書いてあるのがarrに相当する。そんな感じで、表記上は、Linでも全く同じことになって別にSuper使うことに特にメリットはない気がした。実際に、Linで同じことをやってみた

module QComp where
import Complex

class Eq a => Basis a where
 basis :: [a]

instance Basis Bool where
 basis = [True , False]
instance (Basis a,Basis b)=>Basis(a,b) where
 basis = [ (a,b) |a<-basis , b<-basis]
instance (Basis a,Basis b,Basis c) =>Basis(a,b,c) where
 basis = [(a,b,c) | a<-basis,b<-basis,c<-basis]

type Vec a = a -> Complex Double
type Lin a b = a -> Vec b

return :: (Basis a) => a->Vec a
return a =(\b->if (a==b) then 1.0 else 0.0)

bind :: (Basis a) => Vec a -> (a -> Vec b) -> Vec b
bind va f =(\b -> sum[(va a)* ((f a) $ b) |a<-basis])

vplus :: Vec a -> Vec a -> Vec a
vplus v1 v2 = (\a -> (v1 a)+(v2 a))

vminus :: Vec a -> Vec a -> Vec a
vminus v1 v2 = (\a -> (v1 a)-(v2 a))

--tensor product
(<*>) :: Vec a -> Vec b -> Vec(a,b)
v1 <*> v2 =(\(a,b) ->(v1 a) *(v2 b))

--スカラー倍
($*) :: Complex Double -> Vec a -> Vec a
pa $* v = (\a -> pa *(v a))

adjoint :: Lin a b -> Lin b a
adjoint f = (\b -> (\a -> conjugate ((f a) b)))

arr :: (Basis a,Basis b) => (a -> b) -> Lin a b
arr f = (\a -> QComp.return (f a))

(>>>) :: (Basis b,Basis c,Basis d) => Lin b c -> Lin c d -> Lin b d
f >>> g = (\b -> bind (f b) g)

first :: (Basis b,Basis c,Basis d) => Lin b c -> Lin (b,d) (c,d)
first f =(\(b,d) -> (f b) <*> (QComp.return d))

controlled :: Basis a => Lin a a -> Lin (Bool,a) (Bool,a)
controlled f = (\(b1,b2) -> (QComp.return b1) <*> (if b1 then (f b2) else (QComp.return b2)))

--controlled NOT
cnot :: Lin (Bool,Bool) (Bool,Bool)
cnot = controlled (arr not)

--Hadamard transformation
hadamard :: Lin Bool Bool
hadamard False = (1/sqrt 2) $* (vplus (QComp.return True) (QComp.return False))
hadamard True = (1/sqrt 2) $* (vminus (QComp.return True) (QComp.return False))

--なんぞこれ
phase :: Lin Bool Bool
phase False = QComp.return False
phase True = (0 :+ 1) $* (QComp.return True)

cphase = controlled QComp.phase
caphase = controlled (adjoint QComp.phase)

--toffoli Gate(合ってる?
toffoli :: Lin (Bool,Bool,Bool) (Bool,Bool,Bool)
toffoli = 
 arr (\(a0,b0,c0) -> (c0 , (a0,b0))) >>>
 (first hadamard) >>>
 arr (\(c1,(a0,b0)) -> ((b0,c1),a0)) >>>
 (first cphase) >>>
 arr (\((b1,c2),a0) -> ((a0,b1),c2)) >>>
 (first cnot) >>>
 arr (\((a1,b2),c2) -> ((b2,c2),a1)) >>>
 (first caphase) >>>
 arr (\((b3,c3),a1) -> ((a1,b3),c3)) >>>
 (first cnot) >>>
 arr (\((a2,b4),c3) -> ((a2,c3),b4)) >>>
 (first cphase) >>>
 arr (\((a3,c4),b4) -> (c4,(a3,b4))) >>>
 (first hadamard) >>>
 arr (\(c5,(a3,b4)) -> (a3,b4,c5))

で、SuperとLin、どっち使った方がいいのかっていうと、まあ、どっちでもよいと思う。シミュレータとしては、結局やることは行列の計算で、そんでもって、Superにしろ、Linにしろ遅いだろう。なので、結局どっちを使うのかはnotationをどうするかって話に過ぎなくて、そして(型を書かなければ)見掛け上はどっちも変わらない。物理的に、数学的に、どっちの定式化がいいかというと、有限自由度で考えてる限りは、別にどっちでもいいだろう。

Lin a aはユニタリーでない変換を含んでいたり、Superも実際には、Superoperatorでないものを含んでいて、そこらへんは、Haskellでは割とどうしようもない気がする。arrとか使ってコード書くと、実際には実現できない回路が書けてしまう可能性があって、論文の最後の方にごにょごにょ書いてるのは、そういうのをどうするかって話を書いてるんだと思うけど、arrとか使わず、基本ゲートを使って書くようにすれば解決する話で、まあ別段気にするほどの話でもないと思う

あと、Linは、線形代数分かる人にとっては、Arrowを理解するための結構いい例なんじゃないかと思った

7月

■今更

大分前にHaskellで、多項式周りのアルゴリズムを色々実装して遊んでたけど、多項式型は、たぶん、誰でも考えるとおり

data Polynomial a = Const a
                  | Var String
                  | PLUS (Polynomial a) (Polynomial a)
                  | MULT (Polynomial a) (Polynomial a) deriving (Show,Read)

instance (Num a,Ord a) => Num (Polynomial a) where
 p+q = normalize(PLUS p q)
 p*q = normalize(MULT p q)
 p-q = normalize(PLUS p (MULT (Const (negate 1)) q)
 fromInteger n = Const (fromInteger n)

みたいな定義になってる。normalizeは、多項式を適当な標準形に直す関数。

で、がんばってパーサ書いたり、変数に値代入した時は〜とか色々やってたけど、Haskellって

Prelude> (\x y -> 2*x*x+y+3) (Var "x") (Var "y")
PLUS (MULT (Const 2) (MULT (Var "x") (Var "x"))) (PLUS (Var "y") (const 3))
Prelude> (\x y -> 2*x*x+y+3) 4 (Var "hoge")
PLUS (Var "hoge") (Const 35)

とかできてしまうんだよなぁって、気付いた。大したことじゃないけど、私は結構感動した。

逆もできるとよい。つまり、PLUS (MULT (Const 2) (MULT (Var "x") (Var "x"))) (PLUS (Var "y") (const 3))みたいな表現から、(\x y -> 2*x*x+y+3)みたいな式を得たい。式つくるだけなら、

--Language.Haskell.THをimportすること
expandToHs :: Polynomial Integer -> Exp
expandToHs p = LamE (map VarP $ vars p) (mkBody p)
 where
  vars(Const _) = []
  vars(Var x) = [mkName x]
  vars(PLUS p q) = nub $ (vars p)++(vars q)
  vars(MULT p q) = nub $ (vars p)++(vars q)
  mkBody (Const n) = (LitE (IntegerL n))
  mkBody (Var v) = (VarE (mkName v))
  mkBody (PLUS p q) = InfixE (Just (mkBody p)) (VarE (mkName "+")) (Just (mkBody q))
  mkBody (MULT p q) = InfixE (Just (mkBody p)) (VarE (mkName "*")) (Just (mkBody q))

とかでしまい

Prelude> ppr $ expandToHs ((\x y ->2*x*x*x+y*3) (Var "a") (Var "b"))
\a b -> (2 * (a * (a * a))) + (3 * b)

で、evalすればいいんだけど、そうすると、モナドが顔を出すので、なんかイヤ。型に厳しいhaskellでは仕方ないか

■量子こんぷーた

Modeling quantum computing in haskell

http://www.cs.indiana.edu/~sabry/papers/quantum.pdf

別に、こんなん探してたわけではないんだけど、偶然見付けてなんか笑ってしまった。ちなみに読んでない。sabryって、確かFelleisenと、call/ccの公理つくったりしてたような気がする。最近は、量子プログラミング(?)の論文を何本か書いてるようで

http://www.cs.indiana.edu/~sabry/

ぐぐったら、quantum programming languageみたいなキーワードで研究してるひとは結構いるっぽい。Sabryの論文が、量子回路ベースなのに対して、簡単な高級言語+コンパイラ(というか、対応する量子回路を直接生成するっぽいけど)も提唱されている模様。

Compiling a functional quantum programming language

http://www.cs.nott.ac.uk/~gmh/bctcs-slides/grattage.pdf

まだ、まともにマシンもないのに、よくやる。今年頭に、うさんくさいと評判のデモは一応あったけど

http://www.youtube.com/watch?v=VQul2asgXbw

SICPの次の版があるなら、量子コンピュータの章とか設けられたらおもろい

■Superoperatorってなに

#ここに書いたことは、論文チラッと見て妄想で補完しただけなので、全部間違ってるかもしれない

そのSabryの最近の論文で、density matrix+superoperatorで量子力学を定式化すると、superoperatorはArrowの公理満たして、だから、量子回路の記述には、Arrow使うといいんじゃね?みたいなことが書いてあった気がする

superoperatorがArrowの公理満たすってのが何か深い意味がある事実なのかどうかはともかく(単にsuperoperatorの積はsuperoperatorだとか、そういう程度のことを言ってるだけなんじゃないかと思う)、古典的な普通の論理回路でもArrow使うといいよってのは、Hughesが最初から言ってるんで、それはいいことにする。

で、密度行列は分かるけど，superoperatorとか初めて聞いたのだった

http://en.wikipedia.org/wiki/Superoperator

によると、トレースを減らさない、完全正値写像という定義。多分、密度行列の時間発展は、superoperatorになるよってことだと思うんだけど、

・密度行列のトレースって時間発展しても1じゃね？トレースが減少しないって条件は弱いような

・密度行列は正行列だから、時間発展がpositive mapってのは自明だけど、completely positive mapになるって条件は強いような

という2点が腑に落ちない。まぁ、実はsuperoperatorは時間発展じゃないっていうオチかもしれないけど

以下，V,Wを複素ベクトル空間としたとき、Hom(V,V)からHom(W,W)へのsuperoperatorを単に、VからWへのsuperoperatorと呼ぶことにする。論文には、SuperというArrowがでている。例えばSuper (Bool,Bool,Bool) Boolという型は、概念的には、8次元ベクトル空間から、2次元ベクトル空間へのsuperoperatorの全体と考えればよいんだろう。通常の波動関数+ユニタリ変換という定式化に比べてのメリットは、一般にsuperoperatorの定義域と値域(上のVとW)が違っていても構わないってことかな。

問題点。最後の方に、この定義だと、arr(\(x,y)->y)みたいな物理的にありえない計算というか回路も書けてしまう、これ克服するには、線形型備えた「量子プログラミング言語」使うといいんじゃね?というようなことがサラッと書いてあるけど、どうだろうね。

■FPGA最終回

なんかLavaに2種類あるとか書いてある。Xilinx-Lavaってのは、そもそもどこかで入手可能なのかもよく分からんのだけど。まぁいいか

http://www.cs.chalmers.se/~koen/Lava/

それはいいんだけど、とりあえず分かったことは、今のLavaは全然未完成。Lavaのサイトに、思い切り

An implementation of Lava for Xilinx's FPGAs will be made available at the Lava download page shortly. Until there here is a snapshot: lava-0.1b.tar.gz

とか書いてあった

気付くのおせえよ。まあ、現時点でも組合せ回路はがんばれば書けると思うけど，順序回路は今のLavaでは多分書けない気がする。というわけで、Lavaが完成するまで待とう。いつ完成するのか知らないけど。VHDL覚えるとか、自分でLava補完するとかいう選択肢もあるけど、まあめんどくさいな

・variable-freeな全加算器

ふるあだー程度でも、基本関数とコンビネータの組合せで書くのは難しいかもしれないと書いたけど，そんなこともなかった。まぁ、なにを基本関数とするのか、決めてないんでアレだけども

halfAdder = fork2 >-> (and2 <|> xor2)
fullAdder = (halfAdder <|> return) >-> flip' >-> (return <|> halfAdder) >-> flip' >-> (return <|> or2)
 where
  flip' = (flip2 <|> return) >-> switch

flip2 :: (Combinational m bit) => (bit,bit) -> m (bit,bit)
flip2(x,y) = return (y,x)

switch :: (Combinational m bit) => ((bit,bit),bit) -> m (bit , (bit,bit))
switch((x,y),z)=(x,(y,z))

とかで。fork2はLavaに含まれてる関数で、

fork2 a = return (a,a)

という定義。

まあ、flip'が、なんとなくエレガントにいかないというか、こういうのを差し込んでいくと無駄に記述が長くなる。それに

fullAdder ((x,y),z) 
 = do (s1,c1) <- halfAdder (x,y)
      (s2,c2) <- halfAdder (s1,z)
      c3 <- or2(c1,c2)
      return (s2,c3)

の方が分かり易いことは否めない

久しぶりに見て、レイアウトおかしいなぁと思ったら、OS変わったせいだった

■最近学習したもの

Mばっか

・MD(分子動力学)

元の理論というかモデルというか、エネルギー関数が近似というか、経験則でしかないのに、それをがっつり精度あげてみたり、長時間計算したりすることに意味はあるんだろうか

・MO(分子軌道)計算

時間かかりすぎ！

・Mass Spectrometryのデータ処理のアルゴリズム

多分、今後の人生で一切役に立たないであろう知識群

・MPI

まぁ、効率出すのは大変だなぁ。

・Erlang

なんか流行ってたから。なかなかsyntaxに慣れなかったけど、semanticsの論文とか読んでるとおもろい。多分、これ使ってなにか作ることはないだろうけど

・労働基準法

失業したので。世知辛い世の中ですね

■Automated Empirical Optimizations of Software and the ATLAS project

http://citeseer.ist.psu.edu/whaley00automated.html

以前、FFTWをちょろっと眺めたので、ATLASについても眺めて見た。論文は、どうでもいい話が長いので、

http://www.ueda.info.waseda.ac.jp/~ishizaki/pukiwiki/pukiwiki.php?%BF%F4%C3%CD%B7%D7%BB%BB%A5%E9%A5%A4%A5%D6%A5%E9%A5%EAATLAS#content_1_1

コンセプト的には、ATLASは、FFTWの線形代数版みたいな感じで、色んなアーキテクチャに対して、ライブラリを人間が手でチューニングすんのは面倒、自動的にチューニングしようぜ、というような話。けどまあ、FFTWとは、大分やり方が違う

FFTWの方針は、

・原則として、コードレベルでは、"portable"で、マシンレベルの最適化は、plannerが実行時に最適なcodelet(固定したサイズのDFTを計算するコードみたいなもん)の組み合わせを見つけることで行われる

・codelet部分は、genfftが自動生成。codeletは、とにかく、単純に足し算と掛け算の数を減らす

みたいな感じだと思う。

FFTWで、codeletの演算回数を絞るのが、どの程度高速化に寄与してるのかは知らんけど、ATLASの方は、演算量減らすとかはどうでもよくて、メモリアクセスが重要という考えで、うまくレジスタやらキャッシュを使おうというような(FFTWが、そういうのを無視してるというのは言いすぎだけど、でもまぁ直接的には何もしてないし)。そういう最適化は、多分すごい昔から知られてるので、ATLASがやったことってのは、(インストール時に)色んなコードを自動生成して、自動でベンチマークかけて、一番効率のよいのを選ぶ、という人間が手動でやってたことを、ほとんどそのまま自動化した感じっぽくて、あんまり面白みはないような。

あとまあ、ソフトウェアをハードウェアに合わせてチューニングするんじゃなくて、ハード作ろう的な流れはあると思うんだけど。BLAS API(まあ、大体は、行列の掛け算なDGEMM)をFPGAで実装したとかいう論文は、結構あるけど、今のとこ、そういうのが特に利用されてるようでもないし、実際んとこ、汎用CPU+ATLASの組合せと比較して、そんな速いわけでもないか、むしろ遅いというのが実態らしい。遅い理由は、実装が悪いとか、データ転送がボトルネックになってとか色々あるだろうけど。FFTの方は、DSPとかはあるけど、HPCでFFT on FPGAが利用されてるとかって話もあんま聞かないし、どうなんだろう。

■そんなわけで、FPGAでも勉強しよう

<言い訳>まぁ、色々分からないことだらけなのでアレだけども</言い訳>。

方針

・ハードウェア記述言語というと、VHDLとかVerilogとかあるらしいけど、Haskellに埋め込まれたLavaというHDLがあるらしいので、それを使ってみる。理由は特にない

・経済的事情により、当面FPGAボードなし、シミュレータでがんばる

・飽きたらやめる

(たぶん)Lavaの最初の論文

http://citeseer.ist.psu.edu/69503.html

サイト

http://raintown.org/lava/

インストールはmake一発。もしくは、

runhaskell Setup.hs configure
runhaskell Setup.hs build
runhaskell Setup.hs install

GHC6.4.2で何の問題もなし。

論文の頃とは、微妙に色々変わってるようで、結構はまる。まぁ、10年近く前の論文だし。10年経っても、まだバージョン0.1というのは、なんか色々事情があったのかもしれないし、単に中の人が飽きたのかもしれない

あんまり意味ないけど、NANDな例。

import Lava

nand :: Combinational m bit => (bit,bit) -> m bit
nand = and2 >-> inv

nand_top :: Out ()
nand_top
 = do a <- input_bit "a"
      b <- input_bit "b"
      c <- nand (a,b)
      output_bit "c" c 

main
  = do nl <- netlist nand_top
       putXST_VHDL "nandGate" nl

これをコンパイル&実行すると、nandGate.vhdというVHDLファイルが作られる。それが正しいのかどうかは，いまのところ私には分からない(ダメ。まぁ、いくらなんでも、これくらいなら間違えないだろうと信じてる。本体は、

nand=and2>->inv

の部分で、Combinationalとかいう謎のモナドは何かよく分からんがどうでもいいこととする。論文では、Circuitとかいうモナドが登場してたけど、なくなったっぽい(というか、CombinationalとSequentialに分離したのか?)。(>->)は、

(>->) :: (Monad m) => (a -> m b) -> (b -> m c) -> a -> m c

という型。型を見ただけで、何するか大体分かるけど，

f>->g = (\x -> (f x)>>=g)

みたいな定義だと思う(確認してないけど)。Arrowの(>>>)に相当する演算子だと思えばよいっぽい。

mainの部分のnetlistというのは、何してるのかよく分からん。その後のは、なんかVHDLに吐き出してるだけ。VHDLと比べて、劇的に記述量が減るということもない気がする。まぁ、この例だけでは、何とも言えないけど。それでも、私には、VHDLよりはLavaの方がsyntax的に馴染むし、あとまあ、(>->)とか(<|>)みたいなコンビネータが使えるのはよい(<|>は今のLavaでは、消えたっぽい。代わりにpar2とかいう名前になってる)。コンビネータが使えて嬉しいのは、別に記述量が減るからいうより、2つの回路の等価性を機械的に判定できる可能性を示唆しているような気がするので。trivialな例だと、

(a <|> b) >-> (c <|> d) = (a >-> c) <|> (b >-> d)

(a >-> b) >-> c = a >-> (b >-> c)

みたいな。まぁ、一般的にそんなことできるなら、とっくに誰かがやってそうな気もするけど。あと、そもそも、ふるあだーと程度でも、基本関数+コンビネータだけで書くのはしんどいかもしれない

んで、VHDL生成してくれても、それだけではなんもできん。開発ツールには、XilinxのISE WebPackを使おうと思った(使えと書いてあるから)。Linuxのサポートは適当らしいけど、とりあえず、普通にインストールできたと思ったら、何かSynthesizeとかで失敗する。同じことをWindowsでやると、なんか警告はでるけど、成功する。よくわからん。そして、シミュレータを使おうと思っても使えないので、そのうちカーズは考えるのをやめた。

次回以降の目標

・FreeHDLというシミュレータがあるらしいので、そっち試す

・もうちょっと複雑な回路をつくる

・わからんって書いた部分を調べる

キングクリムゾン

■不完全性定理

あんまし興味なかったので、ちゃんと勉強したことなかったけど、勉強しようと思った

http://web.sfc.keio.ac.jp/~mukai/2006-mathlogic/incompleteness.pdf

を読んだ。

たった8ページ(本質的な部分は4ページくらい)で、証明自体に難しいことはないのだけど、

(1)ω-無矛盾の定義がきもい。ω-無矛盾は

自由変数をひとつもつ任意の算術論理式について、全ての自然数nに対してT |- Q(n')ならば、T |- ∀x.Q(x)

という定義でいいのかね。nがメタレベルでの自然数なのに対して、n'は対象言語レベルでの自然数という意味なのだと思うけど、メタレベルでは自然数とか定義されてることになってんの？

(2)背理法とか、ふんだんに使ってるけど、排中律とかは使っていいことになってんの?

とか、そういうあたりで、メタレベルで何を前提としてるのか、よくわからんかった。あと

(3)この証明だと、第一不完全性定理も第二不完全性定理も、PAが無矛盾かつω-無矛盾っていう条件が必要だと思うのだけど(証明可能性が表現可能であるためには、ω-無矛盾でないといけないという以上のことは示してないから)

(4)モーダスポーネンスって一瞬何のことかと思ったよ

上レジュメにある「(第一不完全性定理が)最近、コンピュータを使って証明された」というのは、Boyer-Moore proverによる証明のことかな。

http://www.cs.utexas.edu/users/moore/best-ideas/nqthm/

によると、1986年のことらしいので、最近でもない。この証明では、ペアノ算術じゃなくて、Z2という謎の算術体系に対して証明したらしいのだけど、詳細はよく分からない。最近は、CoqやHOL Lightによる証明もあるらしい

Essential Incompleteness of Arithmetic Verified by Coq

スライド

http://r6.ca/Goedel/TPHOLS2005_incompleteness.pdf

まあ、Coqは多少分かるので、そのうち気が向いたら追って見たい

独立命題について。ゲーデル文とかは、まあ例のための例であって、あんまり面白くないというか、もっと数学者が日常使う体系と日常出会う命題で決定不能な命題がどのくらいあるかという方が興味あるところ。そういうのは、集合論では割と普通にあって、選択公理はZFから独立だし、一般連続体仮説はZFCから独立だし、まあ集合論とか興味ないんで全然知らんけど色々構成されてそう。PAの独立命題では、Paris-Harringtonの定理というのがあるらしい。これは、2階算術では証明可能らしいけど、そもそも数学的に意味のある命題なのかどうかよくわからん。Kruskalのtree theoremとかは、項書換え系とかの論文でたまに見るけど、PAで決定不能だけど、ZFCで証明できる命題らしい。あとまあ、数論には、超越的な方法でしか示されてない命題というのは山ほどあって、その中には普通にPAでは決定不能な命題というのはありそうな気がする。それから、全然知らなかったんだけど、ZFCで決定不能な命題として、Whiteheadの問題とかいうのがあるらしい。これが何か役に立つ定理なのかは分からんけども、数学的にはごく簡単な命題なのに、ZFCでも証明できないってのは結構面白い。ところで、ユークリッド幾何の公理系の独立命題って何かあるんかね

■[雑談]over Turing

人間の意識は計算不可能云々という話を見て思い出したのだけど、計算可能性の定義として、いわゆるチャーチのテーゼは広く受け入れられてる(というか、それ以外の定義はないと思うけど)。けど、チャーチのテーゼが受け入れられてる理由って、単に、チャーチとかゲーデルとかチューリングとかクリーネが計算可能性を定義しようとしたら、等価な定義にたどり着いたという以上の理由はないと思うんで、ほんとにそれでいいの?もっと何かあるんじゃないの?的な疑問を、基礎論とかやってる(まともな)研究者の方が吐露してるのを何年か前に見た。まあ、それは雑談的な話で、だからどうということはなかったけど、プロでも、そういう疑問もつんだなぁとか感心した記憶がある。web上で見たのは確かなんだけど、探し出せない。

あと、こんな論文

Hierarchy of Discrete-Time Dynamical Systems, a Survey

Turing Machineとかを、離散力学系として捉えなおして(という話は、Conwayのライフゲームとかあるので特に変わった話でもないけど)、Turing Machineより強力な離散力学系があるとか、そういう話。まあ、人間の意識はしらんが、宇宙の時間発展は、計算不可能かもしれないとか、そんなことを思った。

■[haskell]Arrow

圏論ばっかやってると飽きる。まあ、モナドが色々な場面でよいインターフェースだという点については、概ね合意が得られてると思うんだけど、Arrowはどうなんかな〜っていうか、そもそも、Arrowとか、まともに使ってるライブラリが殆どないので、よーわからん。モナドがあれば、いつでもArrowが作れるから、そういう意味では一杯あるとも言えるけど

とりあえず、

http://en.wikibooks.org/wiki/Haskell/Understanding_arrows

によると

(1)LL(1)パーサ

Arrowを採用してるというか、John Hughesが"Generalising Monads To Arrows"で、インスパイアされましたって言ってるだけなんだけど。なんか、モナドインターフェースを採用することによって生じる非効率性を解決するためにほげほげ〜的なことが書いてあるけど、とりあえず、Parsecで不満ないしな。ていうか、不満でるほど、パーサなんてしょっちゅう書かないし。あと、ちゃんとした実装が長い間なかったとか。今はPArrowsとかあるらしい

http://www.cs.helsinki.fi/u/ekarttun/PArrows/

PArrow.ToJavaScriptってのが何なのか気になる

(2)Fudgets

Hughesが上の論文で参考にしたって言ってるもう一つのライブラリ。GUIツールキットって書いてるけど、2002年くらいで放置なので、そのまんまではインストールできなかった

(3)Yampa(Arrowised Functional Reactive Programming)

このへんで何か調べてたけど、何か書いてこうと言いながら、全く放置。

(4)Haskell XML ToolBox

http://www.fh-wedel.de/~si/HXmlToolbox/

XML処理用のライブラリは、一応GHC6.6には標準でくっついてた気がするけど、そこを敢えてArrowで、みたいな。

現在のところ、(2)と(3)は放置状態にあるので、まあ何かするなら(1)か(4)だろうなぁ。でも、ArrowXmlって、ArrowApplyのサブクラスで、ArrowApplyからは、ArrowMonadというモナドが作れてしまうので、別にモナドでいいんちゃう?という気がしないでもない。ただ、"Kleisli Arrow"は、自然にArrowChoiceのインスタンスでもあるので、もし、ArrowXmlに入ってるArrowChoiceの構造が、ArrowXmlのArrowMonadのKleisli Arrowから決まるArrowChoice構造と違うなら、意味はあるのかなという気がする。しかし、それを確認するには、Arrowの関係式が多すぎて面倒くさい

■[F#]

最近、やむをえない事情によって、C#を覚えたのだけど、ついでにF#とか使ってみようと思った。ベースになってるのは、Ocamlで、Ocamlの標準ライブラリの大部分を利用できて、.NET Frameworkも簡単に利用できるとか。どうでもいいんだけど、MSには、Erik MeijerとSimon Peyton Jonesがいるんだから、Haskellをベースにしてくれればよかったのに

http://research.microsoft.com/fsharp/release.aspx

v1.1.13.8はなんかインストール時にエラーがでるので、v1.1.12.5をいれた。コンパイラ(fsc)だけでなくて、インタープリタ(fsi)もちゃんとある。Windows以外でも、Monoがあれば動くらしいけど未確認

まあ、Ocamlとか全然分からんのだけど、Ocamlの全ての機能が使えるというわけでもないらしい。けど、

> 2.1 + 3.2;;
val it : float = 5.3
> 2+3;;
val it : int = 5
> 2.1 +. 3.4;;
val it : float = 5.5

浮動小数点演算でも、普通に+が使える。Ocamlだと、+.とかじゃないと型エラーになってダメだった気が。この一点がどうにも嫌いだったので、F#は素晴らしいと思った。まあ、.NETと合わない機能が捨てられてるのかもしれない。そして、Haskellがベースにならなかったのは、あまりにも合わなすぎたとか、そういう事情なんかも

とりあえず、Hello,Worldは

C:\f#>cat hello.fs
print_string "Hello, world!\n";;
printf "Hello, world!\n";;
System.Console.WriteLine "Hello, world!";;
C:\f#>fsc hello.fs
C:\f#>hello
Hello, world!
Hello, world!
Hello, world!

こんな感じ。

Formを開くだけ

#light
open System
open System.Windows.Forms

let form = new Form()
form.Text <- "Hello World Form"
form.Visible <- true
form.Size <- new Drawing.Size(512, 512)

#if COMPILED
[<STAThread>]    
do Application.Run(form)
#endif

#lightとか付けないと、syntax errorになるのは何だろう。caml light構文でパースしてとか、そういう宣言なのかね。全体的にまだ構文がよく分からん。C#使ってて最もあれだったのは、lambdaが欲しいな〜ということだったので(C#3.0には、ラムダ式と呼ばれる物は一応あるけど)、それだけで、F# > C#かな〜と思った

Reflection。foreachは推奨されてないとか警告出るけど、iterだと型推論エラーになるっぽいので。というか、OcamlのArrayと.NETの配列を一緒くたにしてるのがまずいのか。そんくらい面倒見てくれよ〜とか思うけど

#light
open System
open System.Reflection

let x = 3
let members = x.GetType().GetMembers()
do Array.foreach members (fun m->print_endline m.Name)

F#の説明やサンプルは、断片的には探すと結構あるのだけど、どっかに、まとまった資料があるとよい