Paraisoで始める数値流体力学
Paraisoの論文がAcceptしました。確かにこれで、偏微分方程式の解法を記述する数式を処理し、並列仮想マシンを経てCPU/GPU向けのコードを生成し、さらに遺伝的アルゴリズムで自動チューンするまでを一つのシステムとして実証するという一歩を進めることができました。が、現在のParaisoは生後まもない言語で、いろいろとアドホックな箇所が残っています。
インターフェイスも、生成されるコードの速度もまだまだ洗練の余地があり、たくさんのドキュメントを整備していくことで、おそらく誰にも(私にも!)使えない言語という状態を脱してゆきたいと思います。
そのために、数値流体力学の入門記事をParaisoで書いてみようと思います。無論数値流体力学のずっと良い教科書は沢山ありますがParaisoで書いてあるのは初めてなんじゃないかな!インターフェイスというのはみんなで作るものなので、どうかコメントください。
この記事を試すには、ひとまずGitとHaskell Platformが必要です。
また、GPU版を試すにはCUDAが必要になるでしょう。
新しいブランチ
論文の時点からの改良をいれたdonutブランチで、Paraiso 0.3.0.0を作ってゆきます。詳しくはHaddockやParaisoのサイトを参照してください。donutブランチをゲットするにはバージョン管理ソフトウェアgitをインストールした上で次のように入力します。
$ git clone -b donut https://github.com/nushio3/Paraiso.git Cloning into 'Paraiso'... remote: Counting objects: 4743, done. remote: Compressing objects: 100% (1815/1815), done. remote: Total 4743 (delta 2885), reused 4637 (delta 2779) Receiving objects: 100% (4743/4743), 2.65 MiB | 411 KiB/s, done. Resolving deltas: 100% (2885/2885), done. $ cd Paraiso/ $ git branch * donut $
論文執筆時点とくらべ、すでに
- 不評だったLocalとGlobalという用語を、普通にArrayとScalarと呼ぶことにした
- Paraisoが生成するクラスのメソッド名を変更
などの違いが生まれ始めているのでご注意下さい&手間かけてすいません。とまれ、ここで以下のように唱えれば、Paraisoがインストールされるはずです。
Paraiso$ cabal install (...a few minutes...) Registering Paraiso-0.3.0.0... Paraiso$
ここに最初のサンプルがあります。
Paraiso$ cd examples/HelloWorld/
HelloWorld$ cat Generator.hs
#!/usr/bin/env runhaskell {-# LANGUAGE NoImplicitPrelude, OverloadedStrings #-} {-# OPTIONS -Wall #-} import Data.Dynamic import Data.Tensor.TypeLevel import qualified Data.Text.IO as T import Language.Paraiso.Annotation (Annotation) import Language.Paraiso.Generator (generateIO) import qualified Language.Paraiso.Generator.Native as Native import Language.Paraiso.Name import Language.Paraiso.OM import Language.Paraiso.OM.Builder import Language.Paraiso.OM.DynValue as DVal import Language.Paraiso.OM.PrettyPrint (prettyPrintA1) import Language.Paraiso.OM.Realm import qualified Language.Paraiso.OM.Reduce as Reduce import Language.Paraiso.Optimization import NumericPrelude import System.Process (system) -- the names we use table, total, create, tableMaker :: Name table = mkName "table" total = mkName "total" create = mkName "create" tableMaker = mkName "TableMaker" main :: IO () main = do _ <- system "mkdir -p output" T.writeFile "output/OM.txt" $ prettyPrintA1 $ myOM _ <- generateIO mySetup myOM return () mySetup :: Native.Setup Vec2 Int mySetup = (Native.defaultSetup $ Vec :~ 10 :~ 20) { Native.directory = "./dist/" } myOM :: OM Vec2 Int Annotation myOM = optimize O3 $ makeOM tableMaker [] myVars myKernels myVars :: [Named DynValue] myVars = [Named table $ DynValue Array (typeOf (undefined::Int)), Named total $ DynValue Scalar (typeOf (undefined::Int))] myKernels :: [Named (Builder Vec2 Int Annotation ())] myKernels = [Named create createBuilder] createBuilder :: Builder Vec2 Int Annotation () createBuilder = do x <- bind $ loadIndex (undefined::Int) (Axis 0) y <- bind $ loadIndex (undefined::Int) (Axis 1) z <- bind $ x*y store table z store total $ reduce Reduce.Sum z
これは、実行すると九九の表を表示してくれるサンプルです。
HelloWorld$ make runhaskell Generator.hs g++ main.cpp dist/TableMaker.cpp -o main.out HelloWorld$ ./main.out 0 0 0 0 0 0 0 0 0 0 0 1 2 3 4 5 6 7 8 9 0 2 4 6 8 10 12 14 16 18 0 3 6 9 12 15 18 21 24 27 0 4 8 12 16 20 24 28 32 36 0 5 10 15 20 25 30 35 40 45 0 6 12 18 24 30 36 42 48 54 0 7 14 21 28 35 42 49 56 63 0 8 16 24 32 40 48 56 64 72 0 9 18 27 36 45 54 63 72 81 0 10 20 30 40 50 60 70 80 90 0 11 22 33 44 55 66 77 88 99 0 12 24 36 48 60 72 84 96 108 0 13 26 39 52 65 78 91 104 117 0 14 28 42 56 70 84 98 112 126 0 15 30 45 60 75 90 105 120 135 0 16 32 48 64 80 96 112 128 144 0 17 34 51 68 85 102 119 136 153 0 18 36 54 72 90 108 126 144 162 0 19 38 57 76 95 114 133 152 171 total:8550 HelloWorld$