Paraisoで始める数値流体力学

Paraisoの論文がAcceptしました。確かにこれで、偏微分方程式の解法を記述する数式を処理し、並列仮想マシンを経てCPU/GPU向けのコードを生成し、さらに遺伝的アルゴリズムで自動チューンするまでを一つのシステムとして実証するという一歩を進めることができました。が、現在のParaisoは生後まもない言語で、いろいろとアドホックな箇所が残っています。

インターフェイスも、生成されるコードの速度もまだまだ洗練の余地があり、たくさんのドキュメントを整備していくことで、おそらく誰にも(私にも!)使えない言語という状態を脱してゆきたいと思います。

そのために、数値流体力学の入門記事をParaisoで書いてみようと思います。無論数値流体力学のずっと良い教科書は沢山ありますがParaisoで書いてあるのは初めてなんじゃないかな!インターフェイスというのはみんなで作るものなので、どうかコメントください。

この記事を試すには、ひとまずGitHaskell Platformが必要です。
また、GPU版を試すにはCUDAが必要になるでしょう。

新しいブランチ

論文の時点からの改良をいれたdonutブランチで、Paraiso 0.3.0.0を作ってゆきます。詳しくはHaddockParaisoのサイトを参照してください。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$