
\section{\label{formLibrary}The ``form'' library}

\begin{verbatim}
module Lattice (Cell (cellEnum, cellArray,
                      cellAccum, dimension),
                Lattice (MkLattice),
                Form,
                cells, boundary,
                form, accumForm, lattice, rank,
                (!#), (//#),
                listForm, fassocs, findices, felems,
                fmap, fmapi, fzipWith) where

-- The underlying technology of lattices and their
-- forms...  This parallels the `PreludeArray'
-- section for arrays.

infixl 9 !#
infixl 9 //#

-- The `Cell' class is the analogue of the `Ix'
-- class in `PreludeCore'.

class (Eq a) => Cell a where
  cellEnum  :: Int -> [a] -> [a]
  cellArray :: Int -> [a] -> [Assoc a b] -> a -> b
  cellAccum :: Int -> [a] -> [Assoc a b] -> a -> [b]
  dimension :: a -> Int  

-- First field in `Lattice' defines extent of the
-- complex, where required. Second field is the
-- `boundary' function for the complex.

data (Cell a) =>
  Lattice a =
    MkLattice (Int -> [a]) (a -> [Assoc a Int])
data (Cell a) =>
  Form a b  = MkForm (Lattice a) Int (a -> b)

cells     :: (Cell a) => Lattice a -> Int -> [a]
boundary  :: (Cell a, Num b) =>
             Lattice a -> a -> [Assoc a b]


form      :: (Cell a) =>
   Lattice a -> Int -> [Assoc a b] -> Form a b
accumForm :: (Cell a) =>
   Lattice a -> Int -> [Assoc a b] -> Form a [b]
lattice   :: (Cell a) => Form a b -> Lattice a
rank      :: (Cell a) => Form a b -> Int
(!#)      :: (Cell a) => Form a b -> a -> b

(//#)     :: (Cell a) =>
             Form a b -> [Assoc a b] -> Form a b

listForm  :: (Cell a) =>
             (Lattice a) -> Int -> [b] -> Form a b

fassocs   :: (Cell a) => Form a b -> [Assoc a b]
findices  :: (Cell a) => Form a b -> [a]
felems    :: (Cell a) => Form a b -> [b]

fmap      :: (Cell a) => (b -> c) ->
                         Form a b -> Form a c
fmapi     :: (Cell a) => (a -> b -> c) ->
                         Form a b -> Form a c
fzipWith  :: (Cell a) => (b -> c -> d) -> Form a b ->
                         Form a c -> Form a d

-- `Lattice' primitives...

cells    (MkLattice cls _) r = cellEnum r (cls r)
boundary (MkLattice _ bdy) i =
  [j := fromIntegral s | j := s <- bdy i]
 
-- `Form' primitives...

form      l@(MkLattice cls _) r ivs =
    MkForm l r (cellArray r (cls r) ivs)
accumForm l@(MkLattice cls _) r ivs =
    MkForm l r (cellAccum r (cls r) ivs)
lattice (MkForm l _ _) = l
rank    (MkForm _ r _) = r
(!#)    (MkForm _ _ d) = d

-- Higher level operations on Forms...

a //# us   = form (lattice a) (rank a)
    ([i := a !# i
	| i <- findices a \\ [i | i := _ <- us]]
      ++ us)

listForm l r vs =
    form l r (zipWith (:=) (cells l r) vs)

fassocs  a =
    [i := a !# i | i <- cells (lattice a) (rank a)]
findices a = cells (lattice a) (rank a)
felems   a = [a !# i | i <- findices a]

fmap   f a =
    form l r [i := f (a !# i) | i <- cells l r]
    where l = lattice a
          r = rank    a

fmapi  f a =
    form l r [i := f i (a !# i) | i <- cells l r]
    where l = lattice a
          r = rank    a

fzipWith f a b =
    form l r [i := f (a !# i) (b !# i)
                | i <- cells l r]
    where l = lattice a
          r = rank    a
\end{verbatim}

\section{\label{cellInstance}An example instance of the ``Cell'' class}

\begin{verbatim}
module Cell2d (Cell2d (Cell2d0, Cell2d1, Cell2d2)) where

-- A `Cell' instance suitable for representing
-- cells in a 2d grid.

import Lattice (Cell (cellEnum, cellArray, cellAccum))

data Cell2d = Cell2d0 (Int, Int)     |
              Cell2d1 (Int, Int) Int |
              Cell2d2 (Int, Int)       deriving (Text)

ix0 (Cell2d0 i)        = i
ix1 (Cell2d1 (x, y) d) = (x, y, d)
ix2 (Cell2d2 i)        = i

bd0 [c1, c2]    = (ix0 c1, ix0 c2)
bd1 [c1, c2]    = (ix1 c1, ix1 c2)
bd2 [c1, c2]    = (ix2 c1, ix2 c2)

cell0 i         = Cell2d0 i
cell1 (x, y, d) = Cell2d1 (x,y) d
cell2 i         = Cell2d2 i

instance Cell Cell2d where

  cellEnum r cls =
    case r of
      0 -> [cell0 i | i <- range (bd0 cls)]
      1 -> [cell1 i | i <- range (bd1 cls)]
      2 -> [cell2 i | i <- range (bd2 cls)]

  cellArray r cls ivs =
    case r of
      0 -> (array (bd0 cls)
            [ix0 i := v | i := v <- ivs] !) . ix0
      1 -> (array (bd1 cls)
            [ix1 i := v | i := v <- ivs] !) . ix1
      2 -> (array (bd2 cls)
            [ix2 i := v | i := v <- ivs] !) . ix2

  cellAccum r cls ivs =
    case r of
      0 -> (accumArray (flip (:)) [] (bd0 cls)
            [ix0 i := v | i := v <- ivs] !) . ix0
      1 -> (accumArray (flip (:)) [] (bd1 cls)
            [ix1 i := v | i := v <- ivs] !) . ix1
      2 -> (accumArray (flip (:)) [] (bd2 cls)
            [ix2 i := v | i := v <- ivs] !) . ix2
\end{verbatim}

\section{\label{topologyExample}Example definition of lattice geometry}

\begin{verbatim}
module Torus2d (torus2d) where

-- `torus2d' returns an object of type `Lattice Cell2d'
-- implementing a 2 dimensional grid with
-- periodic boundary conditions.

import Lattice (Cell (cellEnum, cellArray, cellAccum),
                Lattice (MkLattice))

import Cell2d (Cell2d (Cell2d0, Cell2d1, Cell2d2))

torus2d :: Int -> Int -> Lattice Cell2d
torus2d lx ly = MkLattice (cls lx ly) (bdy lx ly)

-- Lattice bounds, specified by smallest and largest
-- cell of each dimension:

cls lx ly 0 = [Cell2d0 (0, 0),
               Cell2d0 (lx - 1, ly - 1)]
cls lx ly 1 = [Cell2d1 (0, 0) 0,
               Cell2d1 (lx - 1, ly - 1) 1]
cls lx ly 2 = [Cell2d2 (0, 0),
               Cell2d2 (lx - 1, ly - 1)]

-- Lattice topology (cell boundary):

bdy _ _ (Cell2d0 _) = []

bdy lx ly (Cell2d1 (x, y) d) =
  [Cell2d0 (x, y) :=  -1,
   Cell2d0 (case d of
              0 -> ((x + 1) `mod` lx, y)
              1 -> (x, (y + 1) `mod` ly)) := 1]

bdy lx ly (Cell2d2 (x, y)) =
  [Cell2d1 (x, y) 0 :=  1,
   Cell2d1 ((x + 1) `mod` lx, y) 1  :=  1,
   Cell2d1 (x, (y + 1) `mod` ly) 0  := -1,
   Cell2d1 (x, y) 1 := -1]
\end{verbatim}

\section{Running some of the codes}

     In the main text detailed code was given for evolving
the states of various systems, but initialisation code was
omitted.  As an example, this appendix shows how to prepare
an initial state for the XY model simulation of section
\ref{XYimplementation}.  This is not quite trivial, because
the state incorporates an infinite list of pseudorandom
numbers at each site.

     Suppose $seed$ is an integer used as a random number seed, and $n$
is the size of a square lattice.  Also $theta$ is a normalised
temperature parameter.  In a physicist's simulation $theta$ might be
around unity and $n$ might be somewhere in the range 10-100.  Our
unoptimised Haskell program would be limited to the bottom end of this
range.  To create a lattice with periodic boundary conditions we define
\begin{tabbing}
\verb$  $\= \kill
\>\verb$l   = torus2d $$n$\verb$ $$n$
\end{tabbing}
Then a suitable initial state, fed to \verb$evolve$, is
\begin{tabbing}
\verb$  $\= \kill
\>\verb$ini = listForm l 0 [MkStat $$theta$\verb$ r (1.0 :+ 0.0) |$ \\
\>\verb$                    r <- ranStreams nvars $$seed$\verb$]$ \\
\>\verb$      where nvars = length (felements ini)$
\end{tabbing}
In this initial state all the $u$ variables are set to unity.
\verb$nvars$ is the number of these variables.  For a square lattice it
will just be $n^2$, but for generality it is calculated as the number
of elements of \verb$ini$.  \verb$ranStreams nvars $$seed$ returns a
list of \verb$nvars$ independent random number streams.  Our
leapfrogged linear congruential implementation is given below.  Note
the use of arbitrary precision \verb$Integer$ arithmetic.  Randoms are
scaled to the interval $[0, 1)$.
\begin{verbatim}
module Random (ranStreams) where
 
-- Random number generation.
 
a = 6364136223846793005 :: Integer
b = 1                   :: Integer
r = 2^64                :: Integer
 
next :: (Integer, Integer) -> Integer -> Integer
next (m, c) x = ((m * x) + c) `mod` r
 
intRanStream :: Int -> [Integer]
intRanStream seed = iterate (next (a, b)) (toInteger seed)
 
leap :: Int -> (Integer, Integer) -> (Integer, Integer)
-- Evaluate parameters of n-fold leapfrog generator, given
-- parameters (a, b) of simple linear congruential generator.
leap  0       _     = (1, 0)
leap (n + 1) (a, b) = ((a * m) `mod` r, (b * m + c) `mod` r)
                      where (m, c) = leap n (a, b)
 
intRanStreams :: Int -> Int -> [[Integer]]
intRanStreams n seed = map (iterate (next (leap n (a, b)))) seeds
                       where seeds = take n (intRanStream seed)
 
scale :: (Fractional a) => Integer -> a
scale n = (fromInteger n) / (fromInteger r)
 
ranStreams :: Int -> Int -> [[Float]]
ranStreams n = (map (map scale)) . (intRanStreams n)
\end{verbatim}
