Tuesday, May 22, 2012

Writing your first Orthotope Machine Program

Paraiso is language embedded in Haskell, specialized in the domain of explicit solvers of partial differential equations on uniform mesh. The computation of the domain is represented as the programs of the Orthotope Machine(OM). The Orthotope Machine is design to capture all of the parallelism inherent in the original problem. It is free from the order the elements of the array is stored and accessed. The dependencies between the instructions of the machine are clearly specified and independent instructions can be issued in any order or in parallel. Paraiso is the framework built around this OM; it can translate mathematical descriptions of the PDE-solving algorithm to OM programs; it can convert OM programs to native programs in C++ and CUDA; it can automatically tune them using evolutionary computation. Please refer to Paraiso paper for more details.

But in one short sentence, all you need to do is to create a value of type OM and pass it to the code generator.

Dive into the code

Our first sample code examples/HelloWorld/Generator.hs includes everything you need to tell Paraiso to let him generate code for you. Let's look inside it!

The first few lines are something like #include <stdio.h>.

#!/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 first part of the sample defines a few names for the code we are going to generate. The type Name, defined in the module Language.Paraiso.Name, is just a newtype for Text defined in Data.Text module. We newtype it to distinguish the identifiers from the rest of the text.

-- the names we use
table, total, create, tableMaker :: Name
table = mkName "table"
total = mkName "total"
create = mkName "create"
tableMaker = mkName "TableMaker"

First, let's look at main routine. The first line creates an output folder. The second line is for pretty-printing the data-flow graphs of our first Orthotope Machine, namely myOM. This is just for your interest; it's not a necessary step for the code generation. The third line is the most important; it combines the myOM with mySetup and generates the C++ code.

main :: IO ()
main = do
  _ <- system "mkdir -p output"
  T.writeFile "output/OM.txt" $ prettyPrintA1 $ myOM
  _ <- generateIO mySetup myOM
  return () 

The data type Setup, defined in Language.Paraiso.Generator.Native, is responsible for all kinds of data needed at native-code generation timing but not needed in defining the abstract Orthotope Machine. Here, we specify the size of the actual OM to be generated, and the folder for the generated library. The default target language for code generation is C++.

mySetup :: Native.Setup Vec2 Int
mySetup = 
  (Native.defaultSetup $ Vec :~ 10 :~ 20)
  { Native.directory = "./dist/" 
  }

Next, let's look at the definition of myOM. To define an OM, pass four arguments, namely its name, the list of its global annotations, the list of its static variables, the list of its kernels, to function makeOM.

myOM :: OM Vec2 Int Annotation
myOM = optimize O3 $
  makeOM tableMaker [] myVars myKernels

myOM has two static variables. One is the array for the multiplication table; the other is for calculating the total of the table. They are Array andScalar Realm, respectively.

myVars :: [Named DynValue]
myVars = [Named table $ DynValue Array (typeOf (undefined::Int)),
          Named total $ DynValue Scalar (typeOf (undefined::Int))]

myOM has one kernel, called create. A kernel of an OM is a bunch of computation that is performed at once.

myKernels :: [Named (Builder Vec2 Int Annotation ())]
myKernels = [Named create createBuilder]

Finally, you can make OM kernels by using Builder Monad combinators. Here, let's make a two-dimensional OM kernel that loads x and y indices and multiply them. Let's also calculate the grand-sum of the table.

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

No comments:

Post a Comment