HaskellのaccumArray

Matlabのaccumarrayを良く使っていたので，Haskellにもあるかと思って探してみた．
以下はそのメモ．

Data.Arrayの中にaccumArrayがある
accumArray :: (e->a->e) -> e -> (i,i) -> [(i,a)] -> Array i e
引数1. (e->a->e) は a型 の値をどうやって e型 に蓄積していくかを決める関数
引数2. e は e型の値であり，蓄積の初期値
引数3. (i,i)はインデックスの範囲を決めるタプル
引数4. [(i,a)]は蓄積の対象になるデータ
引数5. Array i e は戻り値．[(i,e)]

ここでの目的からでは1の(e->a->e)は(+)で良い．[(i,a)]のiがインデックスで
ヒストグラムのi番目のビンを表すなら，(i,a)の存在はi番目のビンに相当する
Array i e の i番目のeをe+1という操作で蓄積したいため，aは常に1となるよう
に[(i,a)]を生成すればよい．また，この場合の蓄積の初期値は0で良いので，eは
0を渡しておけばよい．(i,i)はインデックスの範囲を決めるタプルなので，1から
10までのインデックスがあるなら(1,10)となる．すなわち(smallest-index,
largest-index)である．
例）Prelude Data.Array> accumArray (+) 0 (1,7) [(1,1),(3,1),(5,1),(7,1)]
       array (1,7) [(1,1),(2,0),(3,1),(4,0),(5,1),(6,0),(7,1)]

inRange :: Ix a => (a,a) -> a -> Bool も同様にData.Arrayの関数
引数1. (a,a)は範囲を表すタプル．(smallest-value, largest-value)だが，
          valueはIntやCharである必要がある．（Ix型．詳細は:i Ixで）
引数2. a は範囲内に存在するかどうかテストする値を表す
引数3. a が (a,a)内にあればTrueを返し，なければFalseを返す
例）Prelude Data.Array> inRange (1,5) 3
       True
    Prelude Data.Array> inRange (1,5) 7
       False
    Prelude Data.Array> inRange ((1,5),(10,12)) (7,8)
       True

Data.ArrayのArray型はNumeric.LinearAlgebraのfromArray2D
などの関数でMatrix型に変換できる

Haskellで楽に並列処理

Haskellはどうやら計算の並列化に向いるらしい．Matlabのparfor並に楽に並列化できるなら，趣味だけじゃなくて研究でもHaskellで楽しくプログラミングできるかもしれないのでちょっと調べてみた．

そしたらまさにそれらしき情報を発見（http://yunomu.hatenablog.jp/entry/2012/05/12/060238）したので，早速試してみる．

monad-parallelというパッケージを使うらしいのでまずはインストール．cabalを使えば楽勝．

$ cabal update
$ cabal list monad-parallel
　* monad-parallel
    Synopsis: Parallel execution of monadic computations
    Default available version: 0.7.1.1
    Installed versions: [ Not installed ]
    Homepage: http://trac.haskell.org/SCC/wiki/monad-parallel
    License:  GPL

$ cabal install monad-parallel

おわり．Windows（Cygwin）でもUbuntuでも特に問題なし．

使える関数はここを参照→
http://hackage.haskell.org/packages/archive/monad-parallel/0.5.1/doc/html/Control-Monad-Parallel.html

じゃあ次はテストをしてみる．
fibTest.hsという名前で以下を保存．


---------- ここから ----------
fib 0 = 1
fib 1 = 1
fib n = fib (n-2) + fib (n-1)

test n = print $ fib n

main = sequence [(test 40), (test 40), (test 40)]
---------- ここまで ----------


要は適当なフィボナッチ数を計算して表示するって計算を3回やるプログラム．
以下で時間を計ってみる．
$ ghc --make fibTest.hs -O
$ time ./fibTest

次にfibTestPar.hsという名前で以下を保存

---------- ここから ----------
import qualified Control.Monad.Parallel as P
fib 0 = 1
fib 1 = 1
fib n = fib (n-2) + fib (n-1)

test n = print $ fib n

main = P.sequence [(test 40), (test 40), (test 40)]
---------- ここまで ----------


コンパイル時および実行時にオプションが必要なので気を付ける．
$ ghc -threaded --make figTestPar.hs -O
$ time ./fibTestPar +RTS -N

パフォーマンスモニターで見れば並列計算が行われて時間が短くなったのが分かる．
実行時の　+RTS -N　オプションはOSが適当に使うコア数を決めるというオプション．
「2つのコアしか使いたくない」とかがあるなら +RTS -N2 でOK．
（よくわかんないけど，+RTS .... -RTS　の間に並列計算のオプションを挟むってことか？）

他にもいろいろオプションがあるみたいだけど，それは必要になったときに調べます．

hmatrixのインストール

Windows+Cygwin環境でHaskellをやりたくて，前回Haskell PlatformとMeadowのhaskell-modeまではやったけど，hmatrixのインストールで往生してた件，やっと情報を見つけた．英語だけどこれ（https://github.com/AlbertoRuiz/hmatrix/blob/master/INSTALL.md）が良いみたい．

1. gsl-lapack-windows.zipの入手

blas, lapack, gsl-0といったhmatrixに要求されているライブラリが全てまとめてあるとの事．ここ
（https://github.com/downloads/AlbertoRuiz/hmatrix/gsl-lapack-windows.zip）からgsl-lapack-windows.zipをダウンロードできる．以下はダウンロードしたgsl-lapack-windows.zipを/usr/local/srcに置いたものとした説明

2. gsl-lapack-windowsを利用したhmatrixのインストール

$ cd /usr/local/src
$ unzip gsl-lapack-windows.zip
$ cabal update
$ cabal install hmatrix --extra-include-dirs="C:\cygwin\usr\local\src\gsl-lapack-windows" --extra-lib-dirs="C:\cygwin\usr\local\src\gsl-lapack-windows"

今度は「fatal error: gsl/gsl_odeiv2.h: No such file or directory」とか言われる．これはgsl-lapack-windowsがまだ古く，最新のhmatrix-0.14.0.1に含まれる一部のソースファイルがコンパイルできないことが原因．hmatrix-0.13.1.0を使えば良いらしい．

$ cabal install hmatrix-0.13.1.0 --extra-include-dirs="C:\cygwin\usr\local\src\gsl-lapack-windows" --extra-lib-dirs="C:\cygwin\usr\local\src\gsl-lapack-windows"

できたー！

3. hmatrixを使用する前に

喜んだのもつかの間，hmatrixを使ったプログラムをインタプリタで実行しようとするとdllが見つからないと怒られてしまう．gsl-lapack-windowsの中にあるblas.dll，lapack.dll，libgsl-0.dll，libgslcblas-0.dllの４つ全てをロードできる場所にコピーしないといけないのだが，64bitのWindowsの場合は”C:\Windows\System32”の中にコピーしても意味が無い．何とも分かりにくいのだが，64bitのWindowsの場合はdllの32/64bitに関わらず，コピー先は”C:\Windows\SysWOW64”．これは歴史的な経緯から来ているらしいが，とにかく64bit Windowsの場合はdllのコピー先がSysWOW64と覚えておけば良いみたい．

Haskell事初め

研究にしてももプログラミング言語にしても，楽しそうなやつにいろいろ手を出してるから業績も製作物も何ともまとまりが無いことになっちゃってるが，
これはもう性分と諦めて Haskell をやってみようと思う．
（Windows＋Cygwinの環境でやってます）

1. 処理系をインストール

GHCとHugsという処理系が有名らしい．Hugsってやつはメンテが滞っているらしいからGHCを使う事にする．GHC＋お役立ちパッケージをインストールしてくれるHaskell Platformってのをインストールする．
インストーラはここ（http://hackage.haskell.org/platform/index.html）からダウンロードできる．
環境変数も追加してくれるので楽ちん．

2. haskell-modeのインストール

Emacs(Meadow)の設定は非常に重要．ここ（http://www.haskell.org/haskellwiki/Haskell_mode_for_Emacs）を辿って諸々のelispをダウンロードして展開．できたディレクトリをhaskell-modeって名前に変更して~/.emacs.dの中に突っ込む．
あとは，.emacs内で，

(add-to-list 'load-path "~/.emacs.d/haskell-mode")
(require 'haskell-mode)
(require 'haskell-cabal)
(add-to-list 'auto-mode-alist '("\\.hsquot; . haskell-mode))
(add-to-list 'auto-mode-alist '("\\.lhsquot; . literate-haskell-mode))
(add-to-list 'auto-mode-alist '("\\.cabal\\'" . haskell-cabal-mode))

とでもやっておけばとりあえずはOK．その他，いろいろな便利設定はここ（http://d.hatena.ne.jp/kitokitoki/20111217/p1）が詳しいっぽい．

3. Hello World

適当な場所にhello.hsというファイルを作り，中に

main = putStrLn "Hello World!"

という１行を加える．
あとはターミナルで，

$ ghc hello.hs
$ ./hello.exe

でOK．コンパイルせずに実行するためには，

$ ghc -e Main.main hello.hs
あるいは，
$ runghc hello.hs




（追記：hmatrixのインストール）

＊結論は上手くいかないってことなので期待しないように．
ある程度仕事に使えそうじゃないと続かないものが更に続かなくなるので行列演算ライブラリの導入は必須．Haskellにはcabalというシステムがあって簡単にライブラリをインストールできる（現時点ではアンインストールが出来ないので手動でファイルを削除する必要があるらしい）．
早速，

$ cabal update
$ cabal install hmatrix

とやってインストール使用とするものの，blas, lapack, gsl-0がないですよと怒られる．そこでCygwinのsetup.exeを実行し，それらに関する-develのパッケージをインストール．ライブラリとインクルードファイルの場所を指定して再度cabalを実行する．

$ cabal install hmatrix --extra-include-dirs="C:\cygwin\usr\include" --extra-lib-dirs="C:\cygwin\lib"

が，まだgsl-0が無いと言われる．gslらしきものは入れてるので困ってしまうが，しょーがないのでGSLを自分でコンパイルしてインストールしてみる事にする．コンパイルには，autoconf, automake, gcc, libtool, makeが必要なので，持ってない場合はCygwinのsetup.exeでインストール．GSLのソースファイルをここ（http://ftp.jaist.ac.jp/pub/GNU/gsl/）からダウンロードしてきて/usr/local/srcに置く．ダウンロードしたのは最新の1.15．

$ cd /usr/local/src
$ tar xzvf gsl-1.14.tar.gz
$ cd gsl-1.14
$ ./configure
$ make

ここでエラーがでてる．バージョン落としてもダメだし，いろいろ調べてみたがこりゃ難しそう．
Cygwin厄介…なんか萎えるわ…．

Numpyで行列の各行，各列，各要素に関数を適用

PythonでRのapplyみたいな関数が無いか探してたらずいぶん迷ってしまった．
やりたかったのは，numpy.linalg.norm()みたいにベクトルに対して適用する関数を
ある行列の各行ごとに適用することで，Rのapplyならaxisというパラメータをいじれば
適用先が各要素なのか各行なのか各列なのか指定できる．

最初はvectorizeという関数を見つけて，これでできると思ったらaxisみたいなパラメータが無い．
つまり，vectorizeは各要素に関数を適用する場合専用ってこと．
んで，時間を無駄にしてapply_along_axis()をみつけた．

ただそれだけ．

defsetfの中のgensymについて

インヴァージョンを定義するためのdefsetfでハマったのでメモ
今回やりたかったことは，(setf (mref array index1 index2) new-element)とし，accessという名前のマクロで参照されるarrayの要素にnew-elementを束縛すること．
common lispには既にarefという関数があるが，ある行の全ての要素へのアクセスや，行列の部分行列を抜き出しを便利にするためのマクロであるmrefを新たに定義し，それに対するインヴァージョンを定義する．

ここで，mrefは，(mref array rs_re cs_ce)とすれば，arrayのrs+1行目からre行目，cs+1列目からce列目の部分行列を返すもので，(mref array 1 _)とかで1行目の横ベクトルをとってくることもできるように作ってある．

ここで問題になったのはdefsetfの変数補足の問題を減らすために導入されているgensym．
今回の問題の単純な例を示すと，

access-fnにあたるマクロ：
(defmacro aho (x) `(list ',x))
update-fnにあたるマクロ：
(defmacro baka (x y) `(list ',x ',y))

まずは，defsetfの使い方１（短いフォーム）でインヴァージョンを定義すると，
(defsetf aho baka)
で，例えば(setf (aho _) 3)を評価すると(_ 3)となる．これは想定通りの動き．

ところが，defsetfの使い方２（長いフォーム）でインヴァージョンを定義すると，
(defsetf aho (x) (y) `(list ',x ',y))
で，同様に(setf (aho _) 3)を評価すると，「The variable _ is unbound.」というエラーが出る．
defsetfはマクロの定義と同様だという頭でいると，xに与えられた _ というシンボルが評価される訳が無く，`(list ',x ',y)は，(list '_ '3)の様に展開されて，(_ 3)が返ってくるはず！と思いこんでしまう．
しかし，defsetfの中でgensymが使われているという事を思い出して考えてみると，
恐らく，
   (defmacro ..................... (x y)
       (let ((tmp1 (gensym))
              (tmp2 (gensym)))
           `(let* ((,tmp1 ,x)
                     (,tmp2 ,y))
                   ....................
みたいな事になっているので，評価されないと思い込んでたシンボル _ が評価されてしまうことになる．
いやはや，奥が深いねぇ～

Common Lispでのnamed-let

Common Lispにはnamed-letが無い．
そこでマクロの勉強がてらnamed-letを実装
（Schemeみたいにletって名前で統一はできませんでした…）

(defmacro named-let (name bind body)
  `(labels ((,name ,(mapcar #'car bind) ,body))
     (,name ,@(mapcar #'cadr bind))

     ))

アホみたいにバッククオートとアンクオートだけで挑戦してたので時間を無駄に書けてしまった．
「,@」…．ちゃんとググれってことね…

Emacsでソースコードの折りたたみ

hs-minor-mode.elってのがデフォルトで含まれてるらしく，こいつを使えば簡単にソースコードを畳めるようになる．畳む範囲を柔軟に決めれる訳ではないが，ソースコードにタグを打つ必要が無いし，結局これが楽だと思う．

.emacsは以下の通りです．

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;
;;; For folding
;;;

;; C coding style
(add-hook 'c-mode-hook
          '(lambda ()
    (hs-minor-mode 1)))
;; Scheme coding style
(add-hook 'scheme-mode-hook
          '(lambda ()
    (hs-minor-mode 1)))
;; Elisp coding style
(add-hook 'emacs-lisp-mode-hook
          '(lambda ()
    (hs-minor-mode 1)))
;; Lisp coding style
(add-hook 'lisp-mode-hook
          '(lambda ()
    (hs-minor-mode 1)))
;; Python coding style
(add-hook 'python-mode-hook
          '(lambda ()
    (hs-minor-mode 1)))

(define-key
  global-map
  (kbd "C-#") 'hs-toggle-hiding)

SBCL + Emacs(Meadow) + SLIME

Common Lispの設定でどハマりしたのでメモ

1. SBCLのインストール．これについては何も迷う事はない． http://www.sbcl.org/platform-table.htmlからWindows用バイナリを取ってきてC:\SBCL\にインストール（SBCLのWindows版にある不具合については，とりあえず気にしない）
2. Meadowのインストール．上に同じく迷うポイントは無いはず． http://www.meadowy.org/meadow/wiki/%E3%83%80%E3%82%A6%E3%83%B3%E3%83%AD%E3%83%BC%E3%83%89 から現行のsetup-ja.exeをとってきてC:\Meadow\にインストール
3. SLIMEのインストール．これがどハマりポイント．2012/2/21現在のSLIMEにはバグがある模様．ロードバスが通ってるディレクトリ（例えば，~/.emacs.d/）の中でcvs -d :pserver:anonymous:anonymous@common-lisp.net:/project/slime/cvsroot co -D 2009-10-12 slimeを実行して古いヤツをもってくるのが楽．
4. auto-complete.elとac-slime.elのインストール．これはWebを探して適当に持ってきたら良い．
5. .emacsには以下を追加

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;
;;; For Lisp (SBCL)
;;;

(setq inferior-lisp-program "C:/SBCL/sbcl.exe")
(add-to-list 'load-path "C:/cygwin/home/"username"/.emacs.d/slime")
(require 'slime)
(slime-setup)
(slime-setup '(slime-autodoc))
(slime-setup '(slime-repl))
(setq slime-net-coding-system 'utf-8-unix)
(add-hook 'inferior-lisp-mode-hook (lambda () (inferior-slime-mode t)))
(add-hook 'lisp-mode-hook (lambda ()
                            (slime-mode t)
                            (show-paren-mode)))


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;
;;; For auto-complete
;;;
(add-to-list 'load-path "~/.emacs.d/")
(require 'auto-complete-config)
(add-to-list 'ac-dictionary-directories "~/.emacs.d//ac-dict")
(ac-config-default)

;; for slime
(require 'ac-slime)
(add-hook 'slime-mode-hook 'set-up-slime-ac)
(add-hook 'slime-repl-mode-hook 'set-up-slime-ac)
(define-globalized-minor-mode real-global-auto-complete-mode
  auto-complete-mode (lambda ()
         (if (not (minibufferp (current-buffer)))
      (auto-complete-mode 1))))
(real-global-auto-complete-mode t)




こうやって書いてみると何てこと無いね．
これのために数時間無駄にしたと思うとダルくなるわ…