diff --git a/knead.cabal b/knead.cabal
--- a/knead.cabal
+++ b/knead.cabal
@@ -1,12 +1,12 @@
 Name:             knead
-Version:          0.2.3
+Version:          0.3
 License:          BSD3
 License-File:     LICENSE
 Author:           Henning Thielemann <haskell@henning-thielemann.de>
 Maintainer:       Henning Thielemann <haskell@henning-thielemann.de>
 Homepage:         http://hub.darcs.net/thielema/knead/
 Category:         Data Structures
-Synopsis:         Repa array processing using LLVM JIT
+Synopsis:         Repa-like array processing using LLVM JIT
 Description:
   This library processes arrays like @Repa@ and @Accelerate@,
   but it uses the just-in-time compiler of @LLVM@
@@ -15,7 +15,8 @@
   that can be run without a GPU.
   You do not need to care about inlining and strictness annotations,
   because the LLVM code is by default inlined and strict.
-  The package is the basis for an LLVM backend for the @Accelerate@ framework.
+  The package is intended as the basis
+  for an LLVM backend for the @Accelerate@ framework.
   .
   Highlights:
   .
@@ -42,12 +43,13 @@
   .
   The name of the package is inspired by the visualization of typical operations
   like reshaping, collapsing a dimension and extruding another one.
-Tested-With:      GHC==7.4.2, GHC==7.8.4, GHC==8.0.1
+Tested-With:      GHC==7.4.2, GHC==7.8.4
+Tested-With:      GHC==8.4.3, GHC==8.6.1
 Cabal-Version:    >=1.6
 Build-Type:       Simple
 
 Source-Repository this
-  Tag:         0.2.3
+  Tag:         0.3
   Type:        darcs
   Location:    http://hub.darcs.net/thielema/knead/
 
@@ -57,11 +59,12 @@
 
 Library
   Build-Depends:
-    llvm-extra >=0.7.3 && <0.8,
-    llvm-tf >=3.1 && <3.2,
+    llvm-extra >=0.8 && <0.9,
+    llvm-tf >=3.1.1 && <3.2,
     tfp >=1.0 && <1.1,
     storable-tuple >=0.0 && <0.1,
     storable-record >=0.0.3 && <0.1,
+    bool8 >=0.0 && <0.1,
     transformers >=0.3 && <0.6,
     utility-ht >=0.0.1 && <0.1,
     base >=4 && <5
@@ -69,10 +72,11 @@
   GHC-Options:      -Wall
   Hs-Source-Dirs:   src
   Exposed-Modules:
-    Data.Array.Knead.Index.Linear.Int
-    Data.Array.Knead.Index.Linear
-    Data.Array.Knead.Index.Nested.Shape
+    Data.Array.Knead.Shape.Nested
+    Data.Array.Knead.Shape.Cubic
+    Data.Array.Knead.Shape.Cubic.Int
     Data.Array.Knead.Expression
+    Data.Array.Knead.Expression.Vector
     Data.Array.Knead.Parameter
     Data.Array.Knead.Simple.Symbolic
     Data.Array.Knead.Simple.ShapeDependent
diff --git a/src/Data/Array/Knead/Code.hs b/src/Data/Array/Knead/Code.hs
--- a/src/Data/Array/Knead/Code.hs
+++ b/src/Data/Array/Knead/Code.hs
@@ -1,7 +1,7 @@
 {-# LANGUAGE TypeFamilies #-}
 module Data.Array.Knead.Code where
 
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 
 import qualified LLVM.Extra.Multi.Value as MultiValue
 
diff --git a/src/Data/Array/Knead/Expression.hs b/src/Data/Array/Knead/Expression.hs
--- a/src/Data/Array/Knead/Expression.hs
+++ b/src/Data/Array/Knead/Expression.hs
@@ -5,15 +5,15 @@
 import qualified LLVM.Extra.Multi.Value as MultiValue
 import qualified LLVM.Extra.Arithmetic as A
 import qualified LLVM.Extra.Control as C
-import qualified LLVM.Extra.Monad as LMonad
 import qualified LLVM.Core as LLVM
-import LLVM.Extra.Multi.Value (PatternTuple, Decomposed, Atom, atom, )
+import LLVM.Extra.Multi.Value (PatternTuple, Decomposed, Atom, )
 
-import qualified Control.Monad as Monad
+import qualified Control.Monad.HT as Monad
 
 import qualified Data.Tuple.HT as TupleHT
 import qualified Data.Tuple as Tuple
 import Data.Complex (Complex((:+)))
+import Data.Bool8 (Bool8)
 
 import Prelude
    hiding (fst, snd, min, max, zip, unzip, zip3, unzip3,
@@ -47,10 +47,10 @@
 
 instance Value Exp where
    lift0 a = Exp (return a)
-   lift1 f (Exp a) = Exp (Monad.liftM f a)
-   lift2 f (Exp a) (Exp b) = Exp (Monad.liftM2 f a b)
-   lift3 f (Exp a) (Exp b) (Exp c) = Exp (Monad.liftM3 f a b c)
-   lift4 f (Exp a) (Exp b) (Exp c) (Exp d) = Exp (Monad.liftM4 f a b c d)
+   lift1 f (Exp a) = Exp (Monad.lift f a)
+   lift2 f (Exp a) (Exp b) = Exp (Monad.lift2 f a b)
+   lift3 f (Exp a) (Exp b) (Exp c) = Exp (Monad.lift3 f a b c)
+   lift4 f (Exp a) (Exp b) (Exp c) (Exp d) = Exp (Monad.lift4 f a b c d)
 
 
 liftM ::
@@ -65,14 +65,14 @@
     MultiValue.T a -> MultiValue.T b ->
     LLVM.CodeGenFunction r (MultiValue.T c)) ->
    (Exp a -> Exp b -> Exp c)
-liftM2 f (Exp a) (Exp b) = Exp (LMonad.liftR2 f a b)
+liftM2 f (Exp a) (Exp b) = Exp (Monad.liftJoin2 f a b)
 
 liftM3 ::
    (forall r.
     MultiValue.T a -> MultiValue.T b -> MultiValue.T c ->
     LLVM.CodeGenFunction r (MultiValue.T d)) ->
    (Exp a -> Exp b -> Exp c -> Exp d)
-liftM3 f (Exp a) (Exp b) (Exp c) = Exp (LMonad.liftR3 f a b c)
+liftM3 f (Exp a) (Exp b) (Exp c) = Exp (Monad.liftJoin3 f a b c)
 
 
 unliftM1 ::
@@ -94,13 +94,6 @@
 
 
 
-min :: (MultiValue.Real a) => Exp a -> Exp a -> Exp a
-min = liftM2 A.min
-
-max :: (MultiValue.Real a) => Exp a -> Exp a -> Exp a
-max = liftM2 A.max
-
-
 zip :: (Value val) => val a -> val b -> val (a, b)
 zip = lift2 MultiValue.zip
 
@@ -135,10 +128,10 @@
 snd = lift1 MultiValue.snd
 
 mapFst :: (Exp a -> Exp b) -> Exp (a, c) -> Exp (b, c)
-mapFst f = modify (atom, atom) $ TupleHT.mapFst f
+mapFst f = liftM (MultiValue.mapFstF (unliftM1 f))
 
 mapSnd :: (Exp b -> Exp c) -> Exp (a, b) -> Exp (a, c)
-mapSnd f = modify (atom, atom) $ TupleHT.mapSnd f
+mapSnd f = liftM (MultiValue.mapSndF (unliftM1 f))
 
 swap :: (Value val) => val (a, b) -> val (b, a)
 swap = lift1 MultiValue.swap
@@ -160,13 +153,13 @@
 thd3 = lift1 MultiValue.thd3
 
 mapFst3 :: (Exp a0 -> Exp a1) -> Exp (a0,b,c) -> Exp (a1,b,c)
-mapFst3 f = modify (atom, atom, atom) $ TupleHT.mapFst3 f
+mapFst3 f = liftM (MultiValue.mapFst3F (unliftM1 f))
 
 mapSnd3 :: (Exp b0 -> Exp b1) -> Exp (a,b0,c) -> Exp (a,b1,c)
-mapSnd3 f = modify (atom, atom, atom) $ TupleHT.mapSnd3 f
+mapSnd3 f = liftM (MultiValue.mapSnd3F (unliftM1 f))
 
 mapThd3 :: (Exp c0 -> Exp c1) -> Exp (a,b,c0) -> Exp (a,b,c1)
-mapThd3 f = modify (atom, atom, atom) $ TupleHT.mapThd3 f
+mapThd3 f = liftM (MultiValue.mapThd3F (unliftM1 f))
 
 
 modifyMultiValue ::
@@ -375,10 +368,23 @@
 fromRational' = lift0 . MultiValue.fromRational'
 
 
+boolPFrom8 :: Exp Bool8 -> Exp Bool
+boolPFrom8 = lift1 MultiValue.boolPFrom8
+
+bool8FromP :: Exp Bool -> Exp Bool8
+bool8FromP = lift1 MultiValue.bool8FromP
+
+intFromBool8 :: (MultiValue.NativeInteger i ir) => Exp Bool8 -> Exp i
+intFromBool8 = liftM MultiValue.intFromBool8
+
+floatFromBool8 :: (MultiValue.NativeFloating a ar) => Exp Bool8 -> Exp a
+floatFromBool8 = liftM MultiValue.floatFromBool8
+
+
 cmp ::
    (MultiValue.Comparison a) =>
    LLVM.CmpPredicate -> Exp a -> Exp a -> Exp Bool
-cmp ord = liftM2 $ MultiValue.cmp ord
+cmp ord = liftM2 (MultiValue.cmp ord)
 
 infix 4 ==*, /=*, <*, <=*, >*, >=*
 
@@ -392,6 +398,11 @@
 (<=*) = cmp LLVM.CmpLE
 
 
+min, max :: (MultiValue.Real a) => Exp a -> Exp a -> Exp a
+min = liftM2 A.min
+max = liftM2 A.max
+
+
 true, false :: Exp Bool
 true = cons True
 false = cons False
@@ -457,7 +468,8 @@
    abs = liftM MultiValue.abs
    signum = liftM MultiValue.signum
 
-instance (MultiValue.Field a, MultiValue.Real a, MultiValue.RationalConstant a) =>
+instance
+   (MultiValue.Field a, MultiValue.Real a, MultiValue.RationalConstant a) =>
       Fractional (Exp a) where
    fromRational n = lift0 (MultiValue.fromRational' n)
    (/) = liftM2 MultiValue.fdiv
diff --git a/src/Data/Array/Knead/Expression/Vector.hs b/src/Data/Array/Knead/Expression/Vector.hs
new file mode 100644
--- /dev/null
+++ b/src/Data/Array/Knead/Expression/Vector.hs
@@ -0,0 +1,110 @@
+module Data.Array.Knead.Expression.Vector where
+
+import qualified Data.Array.Knead.Expression as Expr
+import Data.Array.Knead.Expression (Exp)
+
+import qualified LLVM.Extra.Multi.Value.Vector as MultiValueVec
+import qualified LLVM.Extra.Multi.Value as MultiValue
+import qualified LLVM.Extra.Multi.Vector as MultiVector
+import qualified LLVM.Core as LLVM
+
+import Prelude hiding (replicate, zip, fst, snd)
+
+
+cons ::
+   (LLVM.Positive n, MultiVector.C a) =>
+   LLVM.Vector n a -> Exp (LLVM.Vector n a)
+cons = Expr.lift0 . MultiValueVec.cons
+
+fst ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b) =>
+   Exp (LLVM.Vector n (a,b)) -> Exp (LLVM.Vector n a)
+fst = Expr.lift1 MultiValueVec.fst
+
+snd ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b) =>
+   Exp (LLVM.Vector n (a,b)) -> Exp (LLVM.Vector n b)
+snd = Expr.lift1 MultiValueVec.snd
+
+swap ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b) =>
+   Exp (LLVM.Vector n (a,b)) -> Exp (LLVM.Vector n (b,a))
+swap = Expr.lift1 MultiValueVec.swap
+
+mapFst ::
+   (Exp (LLVM.Vector n a0) -> Exp (LLVM.Vector n a1)) ->
+   Exp (LLVM.Vector n (a0,b)) -> Exp (LLVM.Vector n (a1,b))
+mapFst f =
+   Expr.liftM
+      (MultiValue.liftM
+         (\(a0,b) -> do
+            MultiValue.Cons a1 <- Expr.unliftM1 f $ MultiValue.Cons a0
+            return (a1,b)))
+
+mapSnd ::
+   (Exp (LLVM.Vector n b0) -> Exp (LLVM.Vector n b1)) ->
+   Exp (LLVM.Vector n (a,b0)) -> Exp (LLVM.Vector n (a,b1))
+mapSnd f =
+   Expr.liftM
+      (MultiValue.liftM
+         (\(a,b0) -> do
+            MultiValue.Cons b1 <- Expr.unliftM1 f $ MultiValue.Cons b0
+            return (a,b1)))
+
+
+fst3 ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b, MultiVector.C c) =>
+   Exp (LLVM.Vector n (a,b,c)) -> Exp (LLVM.Vector n a)
+fst3 = Expr.lift1 MultiValueVec.fst3
+
+snd3 ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b, MultiVector.C c) =>
+   Exp (LLVM.Vector n (a,b,c)) -> Exp (LLVM.Vector n b)
+snd3 = Expr.lift1 MultiValueVec.snd3
+
+thd3 ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b, MultiVector.C c) =>
+   Exp (LLVM.Vector n (a,b,c)) -> Exp (LLVM.Vector n c)
+thd3 = Expr.lift1 MultiValueVec.thd3
+
+
+zip ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b) =>
+   Exp (LLVM.Vector n a) -> Exp (LLVM.Vector n b) ->
+   Exp (LLVM.Vector n (a,b))
+zip = Expr.lift2 MultiValueVec.zip
+
+zip3 ::
+   (LLVM.Positive n, MultiVector.C a, MultiVector.C b, MultiVector.C c) =>
+   Exp (LLVM.Vector n a) -> Exp (LLVM.Vector n b) -> Exp (LLVM.Vector n c) ->
+   Exp (LLVM.Vector n (a,b,c))
+zip3 = Expr.lift3 MultiValueVec.zip3
+
+
+replicate ::
+   (LLVM.Positive n, MultiVector.C a) =>
+   Exp a -> Exp (LLVM.Vector n a)
+replicate = Expr.liftM MultiValueVec.replicate
+
+take ::
+   (LLVM.Positive n, LLVM.Positive m, MultiVector.Select a) =>
+   Exp (LLVM.Vector n a) -> Exp (LLVM.Vector m a)
+take = Expr.liftM MultiValueVec.take
+
+takeRev ::
+   (LLVM.Positive n, LLVM.Positive m, MultiVector.Select a) =>
+   Exp (LLVM.Vector n a) -> Exp (LLVM.Vector m a)
+takeRev = Expr.liftM MultiValueVec.takeRev
+
+
+cmp ::
+   (LLVM.Positive n, MultiVector.Comparison a) =>
+   LLVM.CmpPredicate ->
+   Exp (LLVM.Vector n a) -> Exp (LLVM.Vector n a) -> Exp (LLVM.Vector n Bool)
+cmp ord = Expr.liftM2 (MultiValueVec.cmp ord)
+
+select ::
+   (LLVM.Positive n, MultiVector.Select a) =>
+   Exp (LLVM.Vector n Bool) ->
+   Exp (LLVM.Vector n a) -> Exp (LLVM.Vector n a) -> Exp (LLVM.Vector n a)
+select = Expr.liftM3 MultiValueVec.select
diff --git a/src/Data/Array/Knead/Index/Linear.hs b/src/Data/Array/Knead/Index/Linear.hs
deleted file mode 100644
--- a/src/Data/Array/Knead/Index/Linear.hs
+++ /dev/null
@@ -1,558 +0,0 @@
-{-# LANGUAGE Rank2Types #-}
-{-# LANGUAGE ExistentialQuantification #-}
-{-# LANGUAGE EmptyDataDecls #-}
-{-# LANGUAGE TypeFamilies #-}
-{-# LANGUAGE TypeOperators #-}
-{-# LANGUAGE FlexibleContexts #-}
-{-# LANGUAGE FlexibleInstances #-}
-{-# LANGUAGE UndecidableInstances #-}
-{-# LANGUAGE MultiParamTypeClasses #-}
-module Data.Array.Knead.Index.Linear (
-   C(switch),
-   switchInt,
-   intersect,
-   value,
-   constant,
-   paramWith,
-   tunnel,
-   flattenIndex,
-   peek,
-   poke,
-   computeSize,
-
-   Struct,
-   T(..),
-   Z(Z), z,
-   (:.)((:.)),
-   Shape, shape,
-   Index, index,
-   cons, (#:.),
-   head,
-   tail,
-   switchR,
-   loadMultiValue,
-   storeMultiValue,
-   ) where
-
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
-import qualified Data.Array.Knead.Index.Linear.Int as Index
-
-import qualified Data.Array.Knead.Parameter as Param
-import qualified Data.Array.Knead.Expression as Expr
-import Data.Array.Knead.Expression (Exp, )
-
-import qualified LLVM.Extra.Multi.Value.Memory as MultiValueMemory
-import qualified LLVM.Extra.Multi.Value as MultiValue
-import qualified LLVM.Extra.Arithmetic as A
-import qualified LLVM.Extra.Control as C
-import LLVM.Extra.Multi.Value (Atom, )
-
-import qualified LLVM.Util.Loop as Loop
-import qualified LLVM.Core as LLVM
-
-import qualified Foreign.Storable as St
-import Foreign.Storable.FixedArray (sizeOfArray, )
-import Foreign.Marshal.Array (advancePtr, )
-import Foreign.Ptr (Ptr, castPtr, )
-
-import Control.Monad (liftM2, )
-import Data.Word (Word32, )
-
-import Prelude hiding (min, head, tail, )
-
-
-class C ix where
-   switch ::
-      f Z ->
-      (forall ix0 i. (C ix0, Index.Single i) => f (ix0 :. i)) ->
-      f ix
-
-instance C Z where
-   switch x _ = x
-
-instance (C ix0, Index.Single i) => C (ix0 :. i) where
-   switch _ x = x
-
-
-newtype SwitchInt f ix i = SwitchInt {runSwitchInt :: f (ix :. i)}
-
-switchInt ::
-   (C ix) =>
-   f Z ->
-   (forall ix0. (C ix0) => f (ix0 :. Index.Int)) ->
-   f ix
-switchInt z0 cons0 =
-   switch z0
-      (runSwitchInt $ Index.switchSingle (SwitchInt cons0))
-
-
-newtype Op2 tag sh = Op2 {runOp2 :: Exp (T tag sh) -> Exp (T tag sh) -> Exp (T tag sh)}
-
-intersect :: C sh => Exp (Shape sh) -> Exp (Shape sh) -> Exp (Shape sh)
-intersect =
-   runOp2 $
-   switchInt
-      (Op2 $ \z0 _ -> z0)
-      (Op2 $
-       switchR $ \is i ->
-       switchR $ \js j ->
-          intersect is js #:. Expr.min i j)
-
-
-_value :: (C sh, MultiValue.C sh) => sh -> Exp sh
-_value = Expr.lift0 . MultiValue.cons
-
-
-newtype MakeValue val tag sh = MakeValue {runMakeValue :: T tag sh -> val (T tag sh)}
-
-value :: (C sh, Expr.Value val) => T tag sh -> val (T tag sh)
-value =
-   runMakeValue $
-   switchInt
-      (MakeValue $ \(Cons Z) -> z)
-      (MakeValue $ \(Cons (t:.h)) ->
-         value (Cons t) #:. Expr.lift0 (MultiValue.cons h))
-
-paramWith ::
-   (C sh, Expr.Value val) =>
-   Param.T p (T tag sh) ->
-   (forall parameters.
-    (St.Storable parameters,
-     MultiValueMemory.C parameters) =>
-    (p -> parameters) ->
-    (MultiValue.T parameters -> val (T tag sh)) ->
-    a) ->
-   a
-paramWith p f =
-   case tunnel p of
-      Param.Tunnel get val -> f get (Expr.lift0 . val)
-
-tunnel :: (C sh) => Param.T p (T tag sh) -> Param.Tunnel p (T tag sh)
-tunnel p =
-   case structFieldsPropF p of
-      StructFieldsProp -> Param.tunnel value p
-
-
-data StructFieldsProp sh = LLVM.StructFields (Struct sh) => StructFieldsProp
-
-_structFieldsProp :: (C sh) => f sh -> StructFieldsProp sh
-_structFieldsProp _p = structFieldsRec
-
-structFieldsPropF :: (C sh) => f (g sh) -> StructFieldsProp sh
-structFieldsPropF _p = structFieldsRec
-
-withStructFieldsPropFF ::
-   (C sh) => (StructFieldsProp sh -> f (g (h sh))) -> f (g (h sh))
-withStructFieldsPropFF f = f structFieldsRec
-
-structFieldsRec :: (C sh) => StructFieldsProp sh
-structFieldsRec =
-   switchInt
-      StructFieldsProp
-      (succStructFieldsProp structFieldsRec)
-
-succStructFieldsProp ::
-   StructFieldsProp sh -> StructFieldsProp (sh:.Index.Int)
-succStructFieldsProp StructFieldsProp = StructFieldsProp
-
-
-data Z = Z
-   deriving (Eq, Ord, Read, Show)
-
-
-infixl 3 :., #:.
-
-data tail :. head = !tail :. !head
-   deriving (Eq, Ord, Read, Show)
-
-
-newtype T tag sh = Cons {decons :: sh}
-
-data ShapeTag
-data IndexTag
-
-type Shape = T ShapeTag
-type Index = T IndexTag
-
-shape :: sh -> Shape sh
-shape = Cons
-
-index :: ix -> Index ix
-index = Cons
-
-
-(#:.) :: (Expr.Value val) => val (T tag sh) -> val i -> val (T tag (sh:.i))
-(#:.) = cons
-
-cons :: (Expr.Value val) => val (T tag sh) -> val i -> val (T tag (sh:.i))
-cons =
-   Expr.lift2 $
-      \(MultiValue.Cons t) (MultiValue.Cons h) ->
-         MultiValue.Cons (t,h)
-
-z :: (Expr.Value val) => val (T tag Z)
-z = Expr.lift0 $ MultiValue.Cons ()
-
-head :: (Expr.Value val) => val (T tag (sh:.i)) -> val i
-head = Expr.lift1 $ \(MultiValue.Cons (_t,h)) -> MultiValue.Cons h
-
-tail :: (Expr.Value val) => val (T tag (sh:.i)) -> val (T tag sh)
-tail = Expr.lift1 $ \(MultiValue.Cons (t,_h)) -> MultiValue.Cons t
-
-switchR ::
-   Expr.Value val =>
-   (val (T tag sh) -> val i -> a) -> val (T tag (sh :. i)) -> a
-switchR f ix = f (tail ix) (head ix)
-
-
-instance (tag ~ ShapeTag, sh ~ Z) => Shape.Scalar (T tag sh) where
-   scalar = Expr.lift0 $ MultiValue.Cons ()
-   zeroIndex _ = Expr.lift0 $ MultiValue.Cons ()
-
-
-type family PatternTuple pattern
-type family Decomposed (f :: * -> *) tag pattern
-
-type instance PatternTuple (sh:.s) =
-   PatternTuple sh :. MultiValue.PatternTuple s
-
-type instance Decomposed f tag (sh:.s) =
-   Decomposed f tag sh :. MultiValue.Decomposed f s
-
-type instance PatternTuple (Atom sh) = sh
-
-type instance Decomposed f tag (Atom sh) = f (T tag sh)
-
-
-class
-   (Expr.Composed (Decomposed Exp tag pattern) ~ T tag (PatternTuple pattern)) =>
-      Decompose tag pattern where
-   decompose ::
-      T tag pattern -> Exp (T tag (PatternTuple pattern)) ->
-      Decomposed Exp tag pattern
-
-instance Decompose tag (Atom sh) where
-   decompose (Cons _atom) x = x
-
-instance (Decompose tag sh, Expr.Decompose s) => Decompose tag (sh :. s) where
-   decompose (Cons (psh:.ps)) x =
-      decompose (Cons psh) (tail x) :. Expr.decompose ps (head x)
-
-
-type instance MultiValue.PatternTuple (T tag sh) = T tag (PatternTuple sh)
-
-type instance MultiValue.Decomposed f (T tag sh) = Decomposed f tag sh
-
-
-type family Unwrap sh
-type instance Unwrap (T tag sh) = sh
-
-type family Tag sh
-type instance Tag (T tag sh) = tag
-
-instance
-   (Expr.Compose sh,
-    Expr.Composed sh ~ T (Tag (Expr.Composed sh)) (Unwrap (Expr.Composed sh)),
-    Expr.Compose s) =>
-      Expr.Compose (sh :. s) where
-   type Composed (sh :. s) =
-           T (Tag (Expr.Composed sh))
-             (Unwrap (Expr.Composed sh) :. Expr.Composed s)
-   compose (sh :. s) = cons (Expr.compose sh) (Expr.compose s)
-
-instance (Decompose tag sh) => Expr.Decompose (T tag sh) where
-   decompose = decompose
-
-
-
-instance (C sh) => St.Storable (T tag sh) where
-   sizeOf (Cons sh) = sizeOfArray (rank sh) (0::Word32)
-   alignment (Cons _sh) = St.alignment (0::Word32)
-   poke ptr = poke (castPtr ptr) . decons
-   peek = fmap Cons . peek . castPtr
-
-
-type family Repr (f :: * -> *) sh
-type instance Repr f Z = ()
-type instance Repr f (tail :. head) = (Repr f tail, MultiValue.Repr f head)
-
-instance (C sh) => MultiValue.C (T tag sh) where
-   type Repr f (T tag sh) = Repr f sh
-   cons = value
-   undef = constant $ MultiValue.undef
-   zero = constant $ MultiValue.zero
-   addPhis = addPhis
-   phis = phis
-
-instance (tag ~ ShapeTag, C sh) => Shape.C (T tag sh) where
-   type Index (T tag sh) = Index sh
-   size = fromIntegral . size . decons
-   sizeCode = computeSize
-   intersectCode = Expr.unliftM2 intersect
-   flattenIndexRec sh ix =
-      -- a joint implementation would not be more efficient
-      liftM2 (,)
-         (computeSize sh)
-         (flattenIndex sh ix)
-   loop = loop
-
-
-type family Struct sh
-type instance Struct Z = ()
-type instance Struct (sh :. Index.Int) = (Word32, Struct sh)
-
-instance
-   (C sh, LLVM.StructFields (Struct sh)) =>
-      MultiValueMemory.C (T tag sh) where
-   type Struct (T tag sh) = LLVM.Struct (Struct sh)
-   load = loadMultiValue
-   store = storeMultiValue
-
-loadMultiValue ::
-   (C sh) =>
-   LLVM.Value (Ptr (LLVM.Struct (Struct sh))) ->
-   LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
-loadMultiValue ptr =
-   withStructFieldsPropFF $ \StructFieldsProp ->
-      load =<< castPtrValue ptr
-
-storeMultiValue ::
-   (C sh) =>
-   MultiValue.T (T tag sh) ->
-   LLVM.Value (Ptr (LLVM.Struct (Struct sh))) -> LLVM.CodeGenFunction r ()
-storeMultiValue x ptr =
-   case structFieldsPropF x of
-      StructFieldsProp -> store x =<< castPtrValue ptr
-
-
-newtype FlattenIndex r sh =
-   FlattenIndex {
-      runFlattenIndex ::
-         MultiValue.T (Shape sh) -> MultiValue.T (Index sh) ->
-         LLVM.CodeGenFunction r (LLVM.Value Word32)
-   }
-
-flattenIndex ::
-   (C sh) =>
-   MultiValue.T (Shape sh) -> MultiValue.T (Index sh) ->
-   LLVM.CodeGenFunction r (LLVM.Value Word32)
-flattenIndex =
-   runFlattenIndex $
-   switchInt
-      (FlattenIndex $ \_zerosh _zeroix -> return A.zero)
-      (FlattenIndex $
-         switchR $ \sh (MultiValue.Cons s) ->
-         switchR $ \ix (MultiValue.Cons i) ->
-            A.add i =<< A.mul s =<< flattenIndex sh ix)
-
-
-newtype Rank sh = Rank {runRank :: sh -> Int}
-
-rank :: (C sh) => sh -> Int
-rank =
-   runRank $
-   switch
-      (Rank $ const 0)
-      (Rank $ succ . rank . (\(sh :. _s) -> sh))
-
-
-newtype Peek sh = Peek {runPeek :: Ptr Word32 -> IO sh}
-
-peek :: (C sh) => Ptr Word32 -> IO sh
-peek =
-   runPeek $
-   switchInt
-      (Peek $ const $ return Z)
-      (Peek $ \ptr -> do
-         h <- St.peek ptr
-         t <- peek $ advancePtr ptr 1
-         return (t :. Index.Int h))
-
-
-newtype Poke sh = Poke {runPoke :: Ptr Word32 -> sh -> IO ()}
-
-poke :: (C sh) => Ptr Word32 -> sh -> IO ()
-poke =
-   runPoke $
-   switchInt
-      (Poke $ const $ const $ return ())
-      (Poke $ \ptr (sh :. Index.Int i) -> do
-         St.poke ptr i
-         poke (advancePtr ptr 1) sh)
-
-
-castPtrValue ::
-   (LLVM.StructFields sh) =>
-   LLVM.Value (Ptr (LLVM.Struct sh)) ->
-   LLVM.CodeGenFunction r (LLVM.Value (Ptr Word32))
-castPtrValue = LLVM.bitcast
-
-newtype Load r tag sh =
-   Load {
-      runLoad ::
-         LLVM.Value (Ptr Word32) ->
-         LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
-   }
-
-load ::
-   (C sh) =>
-   LLVM.Value (Ptr Word32) ->
-   LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
-load =
-   runLoad $
-   switchInt
-      (Load $ const $ return z)
-      (Load $ \ptr -> do
-         h <- LLVM.load ptr
-         t <- load =<< A.advanceArrayElementPtr ptr
-         return (t #:. MultiValue.Cons h))
-
-
-newtype Store r tag sh =
-   Store {
-      runStore ::
-         MultiValue.T (T tag sh) ->
-         LLVM.Value (Ptr Word32) ->
-         LLVM.CodeGenFunction r ()
-   }
-
-store ::
-   (C sh) =>
-   MultiValue.T (T tag sh) ->
-   LLVM.Value (Ptr Word32) ->
-   LLVM.CodeGenFunction r ()
-store =
-   runStore $
-   switchInt
-      (Store $ \_z _ptr -> return ())
-      (Store $ switchR $ \sh (MultiValue.Cons k) ptr -> do
-         LLVM.store k ptr
-         store sh =<< A.advanceArrayElementPtr ptr)
-
-
-newtype Size sh =
-   Size {
-      runSize :: sh -> Word32
-   }
-
-size :: (C sh) => sh -> Word32
-size =
-   runSize $
-   switchInt
-      (Size $ \_z -> 1)
-      (Size $ \(sh :. Index.Int k) -> k * size sh)
-
-
-newtype ComputeSize r sh =
-   ComputeSize {
-      runComputeSize ::
-         MultiValue.T (Shape sh) ->
-         LLVM.CodeGenFunction r (LLVM.Value Word32)
-   }
-
-computeSize ::
-   (C sh) =>
-   MultiValue.T (Shape sh) ->
-   LLVM.CodeGenFunction r (LLVM.Value Word32)
-computeSize =
-   runComputeSize $
-   switchInt
-      (ComputeSize $ \_z -> return A.one)
-      (ComputeSize $ switchR $ \sh (MultiValue.Cons k) ->
-         A.mul k =<< computeSize sh)
-
-
-newtype
-   Constant val tag sh =
-      Constant {getConstant :: val Index.Int -> val (T tag sh)}
-
-constant :: (C sh, Expr.Value val) => val Index.Int -> val (T tag sh)
-constant =
-   getConstant $
-   switchInt
-      (Constant $ const z)
-      (Constant $ \x -> constant x #:. x)
-
-
-newtype AddPhis r tag sh =
-   AddPhis {
-      runAddPhis ::
-         LLVM.BasicBlock ->
-         MultiValue.T (T tag sh) ->
-         MultiValue.T (T tag sh) ->
-         LLVM.CodeGenFunction r ()
-   }
-
-addPhis ::
-   (C sh) =>
-   LLVM.BasicBlock ->
-   MultiValue.T (T tag sh) ->
-   MultiValue.T (T tag sh) ->
-   LLVM.CodeGenFunction r ()
-addPhis =
-   runAddPhis $
-   switchInt
-      (AddPhis $ \_ _ _ -> return ())
-      (AddPhis $ \bb ->
-       switchR $ \hx tx ->
-       switchR $ \hy ty ->
-          MultiValue.addPhis bb tx ty >>
-          addPhis bb hx hy)
-
-
-newtype Phis r tag sh =
-   Phis {
-      runPhis ::
-         LLVM.BasicBlock ->
-         MultiValue.T (T tag sh) ->
-         LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
-   }
-
-phis ::
-   (C sh) =>
-   LLVM.BasicBlock ->
-   MultiValue.T (T tag sh) ->
-   LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
-phis =
-   runPhis $
-   switchInt
-      (Phis $ \_ -> return)
-      (Phis $ \bb ->
-       switchR $ \h t ->
-          liftM2 (#:.)
-             (phis bb h)
-             (MultiValue.phis bb t))
-
-
-newtype Loop r state sh =
-   Loop {
-      runLoop ::
-         (MultiValue.T (Index sh) ->
-          state ->
-          LLVM.CodeGenFunction r state) ->
-         MultiValue.T (Shape sh) ->
-         state ->
-         LLVM.CodeGenFunction r state
-   }
-
-loop ::
-   (C sh, Loop.Phi state) =>
-   (MultiValue.T (Index sh) ->
-    state ->
-    LLVM.CodeGenFunction r state) ->
-   MultiValue.T (Shape sh) ->
-   state ->
-   LLVM.CodeGenFunction r state
-loop =
-   runLoop $
-   switchInt
-      (Loop $ \code _z -> code z)
-      (Loop $ \code -> switchR $ \sh (MultiValue.Cons n) ->
-         loop
-            (\ix ptrStart ->
-               fmap fst $
-               C.fixedLengthLoop n (ptrStart, A.zero) $ \(ptr, k) ->
-                  liftM2 (,)
-                     (code (ix #:. MultiValue.Cons k) ptr)
-                     (A.inc k))
-            sh)
diff --git a/src/Data/Array/Knead/Index/Linear/Int.hs b/src/Data/Array/Knead/Index/Linear/Int.hs
deleted file mode 100644
--- a/src/Data/Array/Knead/Index/Linear/Int.hs
+++ /dev/null
@@ -1,59 +0,0 @@
-{-# LANGUAGE TypeFamilies #-}
-module Data.Array.Knead.Index.Linear.Int (
-   Single(..),
-   Int(Int), cons, decons,
-   ) where
-
-import qualified Data.Array.Knead.Expression as Expr
-
-import qualified LLVM.Extra.Multi.Value as MultiValue
-import qualified LLVM.Extra.Arithmetic as A
-
-import Data.Word (Word32, )
-
-import Prelude hiding (Int, head, tail, )
-
-
-newtype Int = Int Word32
-
-cons :: (Expr.Value val) => val Word32 -> val Int
-cons = Expr.lift1 $ \(MultiValue.Cons x) -> MultiValue.Cons x
-
-decons :: (Expr.Value val) => val Int -> val Word32
-decons = Expr.lift1 $ \(MultiValue.Cons x) -> MultiValue.Cons x
-
-
-class Single ix where
-   switchSingle :: f Int -> f ix
-
-instance Single Int where
-   switchSingle x = x
-
-
-instance MultiValue.C Int where
-   type Repr f Int = f Word32
-   cons (Int x) = MultiValue.consPrimitive x
-   undef = MultiValue.undefPrimitive
-   zero = MultiValue.zeroPrimitive
-   phis = MultiValue.phisPrimitive
-   addPhis = MultiValue.addPhisPrimitive
-
-instance MultiValue.Additive Int where
-   add = MultiValue.liftM2 A.add
-   sub = MultiValue.liftM2 A.sub
-   neg = MultiValue.liftM A.neg
-
-instance MultiValue.PseudoRing Int where
-   mul = MultiValue.liftM2 A.mul
-
-instance MultiValue.Real Int where
-   min = MultiValue.liftM2 A.min
-   max = MultiValue.liftM2 A.max
-   abs = MultiValue.liftM A.abs
-   signum = MultiValue.liftM A.signum
-
-instance MultiValue.IntegerConstant Int where
-   fromInteger' = cons . A.fromInteger'
-
-instance MultiValue.Comparison Int where
-   cmp mode = MultiValue.liftM2 $ A.cmp mode
diff --git a/src/Data/Array/Knead/Index/Nested/Shape.hs b/src/Data/Array/Knead/Index/Nested/Shape.hs
deleted file mode 100644
--- a/src/Data/Array/Knead/Index/Nested/Shape.hs
+++ /dev/null
@@ -1,505 +0,0 @@
-{-# LANGUAGE TypeFamilies #-}
-{-# LANGUAGE Rank2Types #-}
-module Data.Array.Knead.Index.Nested.Shape (
-   C(..),
-   value,
-   paramWith,
-   load,
-   intersect,
-   flattenIndex,
-
-   Range(..),
-   Shifted(..),
-
-   Scalar(..),
-   ) where
-
-import qualified Data.Array.Knead.Expression as Expr
-import qualified Data.Array.Knead.Parameter as Param
-import Data.Array.Knead.Expression (Exp, )
-
-import qualified LLVM.Extra.Multi.Value.Memory as MultiMem
-import qualified LLVM.Extra.Multi.Value as MultiValue
-import qualified LLVM.Extra.Arithmetic as A
-import qualified LLVM.Extra.Control as C
-import LLVM.Extra.Multi.Value (atom)
-import LLVM.Extra.Monad (liftR2)
-
-import qualified LLVM.Util.Loop as Loop
-import qualified LLVM.Core as LLVM
-
-import qualified Type.Data.Num.Decimal as TypeNum
-
-import Foreign.Storable
-         (Storable, sizeOf, alignment, poke, peek, pokeElemOff, peekElemOff)
-import Foreign.Ptr (Ptr, castPtr)
-
-import Data.Word (Word8, Word16, Word32, Word64)
-import Data.Int (Int8, Int16, Int32, Int64)
-
-import qualified Control.Monad.HT as Monad
-import Control.Applicative ((<$>))
-
-
-value :: (C sh, Expr.Value val) => sh -> val sh
-value = Expr.lift0 . MultiValue.cons
-
-paramWith ::
-   (Storable b, MultiMem.C b, Expr.Value val) =>
-   Param.T p b ->
-   (forall parameters.
-    (Storable parameters, MultiMem.C parameters) =>
-    (p -> parameters) ->
-    (MultiValue.T parameters -> val b) ->
-    a) ->
-   a
-paramWith p f =
-   Param.withMulti p (\get val -> f get (Expr.lift0 . val))
-
-load ::
-   (MultiMem.C sh) =>
-   f sh -> LLVM.Value (Ptr (MultiMem.Struct sh)) ->
-   LLVM.CodeGenFunction r (MultiValue.T sh)
-load _ = MultiMem.load
-
-intersect :: (C sh) => Exp sh -> Exp sh -> Exp sh
-intersect = Expr.liftM2 intersectCode
-
-flattenIndex ::
-   (C sh) =>
-   MultiValue.T sh -> MultiValue.T (Index sh) ->
-   LLVM.CodeGenFunction r (LLVM.Value Word32)
-flattenIndex sh ix =
-   fmap snd $ flattenIndexRec sh ix
-
-class (MultiValue.C sh) => C sh where
-   type Index sh :: *
-   {-
-   It would be better to restrict zipWith to matching shapes
-   and turn shape intersection into a bound check.
-   -}
-   intersectCode ::
-      MultiValue.T sh -> MultiValue.T sh ->
-      LLVM.CodeGenFunction r (MultiValue.T sh)
-   sizeCode ::
-      MultiValue.T sh ->
-      LLVM.CodeGenFunction r (LLVM.Value Word32)
-   size :: sh -> Int
-   {- |
-   Result is @(size, flattenedIndex)@.
-   @size@ must equal the result of 'sizeCode'.
-   We use this for sharing intermediate results.
-   -}
-   flattenIndexRec ::
-      MultiValue.T sh -> MultiValue.T (Index sh) ->
-      LLVM.CodeGenFunction r (LLVM.Value Word32, LLVM.Value Word32)
-   loop ::
-      (Index sh ~ ix, Loop.Phi state) =>
-      (MultiValue.T ix -> state -> LLVM.CodeGenFunction r state) ->
-      MultiValue.T sh -> state -> LLVM.CodeGenFunction r state
-
-
-instance C () where
-   type Index () = ()
-   intersectCode _ _ = return $ MultiValue.cons ()
-   sizeCode _ = return A.one
-   size _ = 1
-   flattenIndexRec _ _ = return (A.one, A.zero)
-   loop = id
-
-
-class C sh => Scalar sh where
-   scalar :: (Expr.Value val) => val sh
-   zeroIndex :: (Expr.Value val) => f sh -> val (Index sh)
-
-instance Scalar () where
-   scalar = Expr.lift0 $ MultiValue.Cons ()
-   zeroIndex _ = Expr.lift0 $ MultiValue.Cons ()
-
-
-loopPrimitive ::
-   (MultiValue.Repr LLVM.Value j ~ LLVM.Value j,
-    Num j, LLVM.IsConst j, LLVM.IsInteger j,
-    LLVM.CmpRet j, LLVM.CmpResult j ~ Bool,
-    MultiValue.Additive i, MultiValue.IntegerConstant i,
-    Loop.Phi state) =>
-   (MultiValue.T i -> state -> LLVM.CodeGenFunction r state) ->
-   MultiValue.T j -> state -> LLVM.CodeGenFunction r state
-loopPrimitive code (MultiValue.Cons n) ptrStart =
-   loopStart code n MultiValue.zero ptrStart
-
-loopStart ::
-   (Num j, LLVM.IsConst j, LLVM.IsInteger j,
-    LLVM.CmpRet j, LLVM.CmpResult j ~ Bool,
-    MultiValue.Additive i, MultiValue.IntegerConstant i,
-    Loop.Phi state) =>
-   (MultiValue.T i -> state -> LLVM.CodeGenFunction r state) ->
-   LLVM.Value j ->
-   MultiValue.T i -> state -> LLVM.CodeGenFunction r state
-loopStart code n start ptrStart =
-   fmap fst $
-   C.fixedLengthLoop n (ptrStart, start) $ \(ptr, k) ->
-      Monad.lift2 (,) (code k ptr) (MultiValue.inc k)
-
-instance C Word32 where
-   type Index Word32 = Word32
-   intersectCode = MultiValue.min
-   sizeCode (MultiValue.Cons n) = return n
-   size = fromIntegral
-   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) = return (n, i)
-   loop = loopPrimitive
-
-instance C Word64 where
-   type Index Word64 = Word64
-   intersectCode = MultiValue.min
-   sizeCode (MultiValue.Cons n) = LLVM.trunc n
-   size = fromIntegral
-   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
-      Monad.lift2 (,) (LLVM.trunc n) (LLVM.trunc i)
-   loop = loopPrimitive
-
-
-{- |
-Array dimensions and indexes cannot be negative,
-but computations in indices may temporarily yield negative values
-or we want to add negative values to indices.
-
-Maybe we should better have type Index Word64 = Int64?
--}
-unsigned8 :: LLVM.Value Int8 -> LLVM.CodeGenFunction r (LLVM.Value Word8)
-unsigned8 = LLVM.bitcast
-
-instance C Int8 where
-   type Index Int8 = Int8
-   intersectCode = MultiValue.min
-   sizeCode (MultiValue.Cons n) = LLVM.ext =<< unsigned8 n
-   size = fromIntegral
-   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
-      Monad.lift2 (,) (LLVM.ext =<< unsigned8 n) (LLVM.ext =<< unsigned8 i)
-   loop = loopPrimitive
-
-unsigned16 :: LLVM.Value Int16 -> LLVM.CodeGenFunction r (LLVM.Value Word16)
-unsigned16 = LLVM.bitcast
-
-instance C Int16 where
-   type Index Int16 = Int16
-   intersectCode = MultiValue.min
-   sizeCode (MultiValue.Cons n) = LLVM.ext =<< unsigned16 n
-   size = fromIntegral
-   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
-      Monad.lift2 (,) (LLVM.ext =<< unsigned16 n) (LLVM.ext =<< unsigned16 i)
-   loop = loopPrimitive
-
-instance C Int32 where
-   type Index Int32 = Int32
-   intersectCode = MultiValue.min
-   sizeCode (MultiValue.Cons n) = LLVM.bitcast n
-   size = fromIntegral
-   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
-      Monad.lift2 (,) (LLVM.bitcast n) (LLVM.bitcast i)
-   loop = loopPrimitive
-
-instance C Int64 where
-   type Index Int64 = Int64
-   intersectCode = MultiValue.min
-   sizeCode (MultiValue.Cons n) = LLVM.trunc n
-   size = fromIntegral
-   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
-      Monad.lift2 (,) (LLVM.trunc n) (LLVM.trunc i)
-   loop = loopPrimitive
-
-
-{- |
-'Range' denotes an inclusive range like
-those of the Haskell 98 standard @Array@ type from the @array@ package.
-E.g. the shape type @(Range Int32, Range Int64)@
-is equivalent to the ix type @(Int32, Int64)@ for @Array@s.
--}
-data Range n = Range n n
-
-singletonRange :: n -> Range n
-singletonRange n = Range n n
-
-
-{-# INLINE castToElemPtr #-}
-castToElemPtr :: Ptr (f a) -> Ptr a
-castToElemPtr = castPtr
-
--- cf. sample-frame:Stereo
-instance Storable n => Storable (Range n) where
-   {-# INLINE sizeOf #-}
-   {-# INLINE alignment #-}
-   {-# INLINE peek #-}
-   {-# INLINE poke #-}
-   sizeOf ~(Range l r) = sizeOf l + mod (- sizeOf l) (alignment r) + sizeOf r
-   alignment ~(Range l _) = alignment l
-   poke p (Range l r) =
-      let q = castToElemPtr p
-      in  poke q l >> pokeElemOff q 1 r
-   peek p =
-      let q = castToElemPtr p
-      in  Monad.lift2 Range (peek q) (peekElemOff q 1)
-
-
-class
-   (MultiValue.Additive n, MultiValue.Real n, MultiValue.IntegerConstant n) =>
-      ToSize n where
-   toSize :: MultiValue.T n -> LLVM.CodeGenFunction r (LLVM.Value Word32)
-
-instance ToSize Word32 where toSize (MultiValue.Cons n) = LLVM.adapt n
-instance ToSize Word64 where toSize (MultiValue.Cons n) = LLVM.adapt n
-instance ToSize Int32 where toSize (MultiValue.Cons n) = LLVM.bitcast n
-instance ToSize Int64 where toSize (MultiValue.Cons n) = LLVM.trunc n
-
-rangeSize ::
-   (ToSize n) =>
-   Range (MultiValue.T n) -> LLVM.CodeGenFunction r (LLVM.Value Word32)
-rangeSize (Range from to) =
-   toSize =<< MultiValue.inc =<< MultiValue.sub to from
-
-unzipRange :: MultiValue.T (Range n) -> Range (MultiValue.T n)
-unzipRange (MultiValue.Cons (Range from to)) =
-   Range (MultiValue.Cons from) (MultiValue.Cons to)
-
-zipRange :: MultiValue.T n -> MultiValue.T n -> MultiValue.T (Range n)
-zipRange (MultiValue.Cons from) (MultiValue.Cons to) =
-   MultiValue.Cons (Range from to)
-
-instance (MultiValue.C n) => MultiValue.C (Range n) where
-   type Repr f (Range n) = Range (MultiValue.Repr f n)
-   cons (Range from to) = zipRange (MultiValue.cons from) (MultiValue.cons to)
-   undef = MultiValue.compose $ singletonRange MultiValue.undef
-   zero = MultiValue.compose $ singletonRange MultiValue.zero
-   phis bb a =
-      case unzipRange a of
-         Range a0 a1 ->
-            Monad.lift2 zipRange (MultiValue.phis bb a0) (MultiValue.phis bb a1)
-   addPhis bb a b =
-      case (unzipRange a, unzipRange b) of
-         (Range a0 a1, Range b0 b1) ->
-            MultiValue.addPhis bb a0 b0 >>
-            MultiValue.addPhis bb a1 b1
-
-type instance
-   MultiValue.Decomposed f (Range pn) =
-      Range (MultiValue.Decomposed f pn)
-type instance
-   MultiValue.PatternTuple (Range pn) =
-      Range (MultiValue.PatternTuple pn)
-
-instance (MultiValue.Compose n) => MultiValue.Compose (Range n) where
-   type Composed (Range n) = Range (MultiValue.Composed n)
-   compose (Range from to) =
-      zipRange (MultiValue.compose from) (MultiValue.compose to)
-
-instance (MultiValue.Decompose pn) => MultiValue.Decompose (Range pn) where
-   decompose (Range pfrom pto) rng =
-      case unzipRange rng of
-         Range from to ->
-            Range
-               (MultiValue.decompose pfrom from)
-               (MultiValue.decompose pto to)
-
-instance (MultiMem.C n) => MultiMem.C (Range n) where
-   type Struct (Range n) = PairStruct n
-   decompose = fmap (uncurry zipRange) . decomposeGen
-   compose x = case unzipRange x of Range n m -> composeGen n m
-
-instance (Integral n, ToSize n) => C (Range n) where
-   type Index (Range n) = n
-   intersectCode =
-      MultiValue.modifyF2 (singletonRange atom) (singletonRange atom) $
-            \(Range fromN toN) (Range fromM toM) ->
-         Monad.lift2 Range (MultiValue.max fromN fromM) (MultiValue.min toN toM)
-   sizeCode = rangeSize . unzipRange
-   size (Range from to) = fromIntegral $ to-from+1
-   flattenIndexRec rngValue i =
-      case unzipRange rngValue of
-         rng@(Range from _to) ->
-            Monad.lift2 (,) (rangeSize rng) (toSize =<< MultiValue.sub i from)
-   loop code rngValue ptrStart =
-      case unzipRange rngValue of
-         rng@(Range from _to) -> do
-            {-
-            FIXME: rangeSize converts to Word32 which is overly restrictive here.
-            -}
-            n <- rangeSize rng
-            loopStart code n from ptrStart
-
-
-{- |
-'Shifted' denotes a range defined by the start index and the length.
--}
-data Shifted n = Shifted {shiftedOffset, shiftedSize :: n}
-
-singletonShifted :: n -> Shifted n
-singletonShifted n = Shifted n n
-
-
--- cf. sample-frame:Stereo
-instance Storable n => Storable (Shifted n) where
-   {-# INLINE sizeOf #-}
-   {-# INLINE alignment #-}
-   {-# INLINE peek #-}
-   {-# INLINE poke #-}
-   sizeOf ~(Shifted l n) = sizeOf l + mod (- sizeOf l) (alignment n) + sizeOf n
-   alignment ~(Shifted l _) = alignment l
-   poke p (Shifted l n) =
-      let q = castToElemPtr p
-      in  poke q l >> pokeElemOff q 1 n
-   peek p =
-      let q = castToElemPtr p
-      in  Monad.lift2 Shifted (peek q) (peekElemOff q 1)
-
-
-unzipShifted :: MultiValue.T (Shifted n) -> Shifted (MultiValue.T n)
-unzipShifted (MultiValue.Cons (Shifted from to)) =
-   Shifted (MultiValue.Cons from) (MultiValue.Cons to)
-
-zipShifted :: MultiValue.T n -> MultiValue.T n -> MultiValue.T (Shifted n)
-zipShifted (MultiValue.Cons from) (MultiValue.Cons to) =
-   MultiValue.Cons (Shifted from to)
-
-instance (MultiValue.C n) => MultiValue.C (Shifted n) where
-   type Repr f (Shifted n) = Shifted (MultiValue.Repr f n)
-   cons (Shifted offset len) =
-      zipShifted (MultiValue.cons offset) (MultiValue.cons len)
-   undef = MultiValue.compose $ singletonShifted MultiValue.undef
-   zero = MultiValue.compose $ singletonShifted MultiValue.zero
-   phis bb a =
-      case unzipShifted a of
-         Shifted a0 a1 ->
-            Monad.lift2 zipShifted
-               (MultiValue.phis bb a0) (MultiValue.phis bb a1)
-   addPhis bb a b =
-      case (unzipShifted a, unzipShifted b) of
-         (Shifted a0 a1, Shifted b0 b1) ->
-            MultiValue.addPhis bb a0 b0 >>
-            MultiValue.addPhis bb a1 b1
-
-type instance
-   MultiValue.Decomposed f (Shifted pn) =
-      Shifted (MultiValue.Decomposed f pn)
-type instance
-   MultiValue.PatternTuple (Shifted pn) =
-      Shifted (MultiValue.PatternTuple pn)
-
-instance (MultiValue.Compose n) => MultiValue.Compose (Shifted n) where
-   type Composed (Shifted n) = Shifted (MultiValue.Composed n)
-   compose (Shifted offset len) =
-      zipShifted (MultiValue.compose offset) (MultiValue.compose len)
-
-instance (MultiValue.Decompose pn) => MultiValue.Decompose (Shifted pn) where
-   decompose (Shifted poffset plen) rng =
-      case unzipShifted rng of
-         Shifted offset len ->
-            Shifted
-               (MultiValue.decompose poffset offset)
-               (MultiValue.decompose plen len)
-
-instance (MultiMem.C n) => MultiMem.C (Shifted n) where
-   type Struct (Shifted n) = PairStruct n
-   decompose = fmap (uncurry zipShifted) . decomposeGen
-   compose x = case unzipShifted x of Shifted n m -> composeGen n m
-
-type PairStruct n = LLVM.Struct (MultiMem.Struct n, (MultiMem.Struct n, ()))
-
-decomposeGen ::
-   (MultiMem.C n) =>
-   LLVM.Value (PairStruct n) ->
-   LLVM.CodeGenFunction r (MultiValue.T n, MultiValue.T n)
-decomposeGen nm =
-   Monad.lift2 (,)
-      (MultiMem.decompose =<< LLVM.extractvalue nm TypeNum.d0)
-      (MultiMem.decompose =<< LLVM.extractvalue nm TypeNum.d1)
-
-composeGen ::
-   (MultiMem.C n) =>
-   MultiValue.T n -> MultiValue.T n ->
-   LLVM.CodeGenFunction r (LLVM.Value (PairStruct n))
-composeGen n m = do
-   sn <- MultiMem.compose n
-   sm <- MultiMem.compose m
-   rn <- LLVM.insertvalue (LLVM.value LLVM.undef) sn TypeNum.d0
-   LLVM.insertvalue rn sm TypeNum.d1
-
-
-instance (Integral n, ToSize n) => C (Shifted n) where
-   type Index (Shifted n) = n
-   intersectCode =
-      MultiValue.modifyF2 (singletonShifted atom) (singletonShifted atom) $
-            \(Shifted offsetN lenN) (Shifted offsetM lenM) -> do
-         offset <- MultiValue.max offsetN offsetM
-         endN <- MultiValue.add offsetN lenN
-         endM <- MultiValue.add offsetM lenM
-         end <- MultiValue.min endN endM
-         Shifted offset <$> MultiValue.sub end offset
-   sizeCode = toSize . shiftedSize . unzipShifted
-   size (Shifted _offset len) = fromIntegral len
-   flattenIndexRec shapeValue i =
-      case unzipShifted shapeValue of
-         Shifted offset len ->
-            Monad.lift2 (,) (toSize len) (toSize =<< MultiValue.sub i offset)
-   loop code rngValue ptrStart =
-      case unzipShifted rngValue of
-         Shifted from len -> do
-            n <- toSize len
-            loopStart code n from ptrStart
-
-
-
-instance (C n, C m) => C (n,m) where
-   type Index (n,m) = (Index n, Index m)
-   intersectCode a b =
-      case (MultiValue.unzip a, MultiValue.unzip b) of
-         ((an,am), (bn,bm)) ->
-            Monad.lift2 MultiValue.zip
-               (intersectCode an bn)
-               (intersectCode am bm)
-   sizeCode nm =
-      case MultiValue.unzip nm of
-         (n,m) -> liftR2 A.mul (sizeCode n) (sizeCode m)
-   size (n,m) = size n * size m
-   flattenIndexRec nm ij =
-      case (MultiValue.unzip nm, MultiValue.unzip ij) of
-         ((n,m), (i,j)) -> do
-            (ns, il) <- flattenIndexRec n i
-            (ms, jl) <- flattenIndexRec m j
-            Monad.lift2 (,)
-               (A.mul ns ms)
-               (A.add jl =<< A.mul ms il)
-   loop code nm =
-      case MultiValue.unzip nm of
-         (n,m) -> loop (\i -> loop (\j -> code (MultiValue.zip i j)) m) n
-
-instance (C n, C m, C l) => C (n,m,l) where
-   type Index (n,m,l) = (Index n, Index m, Index l)
-   intersectCode a b =
-      case (MultiValue.unzip3 a, MultiValue.unzip3 b) of
-         ((ai,aj,ak), (bi,bj,bk)) ->
-            Monad.lift3 MultiValue.zip3
-               (intersectCode ai bi)
-               (intersectCode aj bj)
-               (intersectCode ak bk)
-   sizeCode nml =
-      case MultiValue.unzip3 nml of
-         (n,m,l) ->
-            liftR2 A.mul (sizeCode n) $
-            liftR2 A.mul (sizeCode m) (sizeCode l)
-   size (n,m,l) = size n * size m * size l
-   flattenIndexRec nml ijk =
-      case (MultiValue.unzip3 nml, MultiValue.unzip3 ijk) of
-         ((n,m,l), (i,j,k)) -> do
-            (ns, il) <- flattenIndexRec n i
-            (ms, jl) <- flattenIndexRec m j
-            x0 <- A.add jl =<< A.mul ms il
-            (ls, kl) <- flattenIndexRec l k
-            x1 <- A.add kl =<< A.mul ls x0
-            sz <- A.mul ns =<< A.mul ms ls
-            return (sz, x1)
-   loop code nml =
-      case MultiValue.unzip3 nml of
-         (n,m,l) ->
-             loop (\i -> loop (\j -> loop (\k ->
-                code (MultiValue.zip3 i j k))
-             l) m) n
diff --git a/src/Data/Array/Knead/Parameter.hs b/src/Data/Array/Knead/Parameter.hs
--- a/src/Data/Array/Knead/Parameter.hs
+++ b/src/Data/Array/Knead/Parameter.hs
@@ -90,6 +90,7 @@
 that use the unit type for constants
 and thus save space in parameter structures.
 -}
+{-# INLINE withTuple #-}
 withTuple ::
    (Storable tuple, Class.MakeValueTuple tuple,
     Class.ValueTuple tuple ~ value, Memory.C value) =>
@@ -105,6 +106,7 @@
 withTuple (Constant a) f = f (const ()) (\() -> Class.valueTupleOf a)
 withTuple (Variable v) f = f v id
 
+{-# INLINE withMulti #-}
 withMulti ::
    (Storable b, MultiValueMemory.C b) =>
    T p b ->
@@ -117,6 +119,7 @@
    a
 withMulti = with MultiValue.cons
 
+{-# INLINE with #-}
 with ::
    (Storable b, MultiValueMemory.C b) =>
    (b -> MultiValue.T b) ->
diff --git a/src/Data/Array/Knead/Parameterized/Physical.hs b/src/Data/Array/Knead/Parameterized/Physical.hs
--- a/src/Data/Array/Knead/Parameterized/Physical.hs
+++ b/src/Data/Array/Knead/Parameterized/Physical.hs
@@ -8,7 +8,7 @@
    the,
    render,
    renderShape,
-   mapAccumL,
+   mapAccumLSimple,
    foldOuterL,
    scatter,
    scatterMaybe,
@@ -20,7 +20,7 @@
 import qualified Data.Array.Knead.Simple.Physical as Phys
 import qualified Data.Array.Knead.Simple.Private as Core
 import qualified Data.Array.Knead.Parameter as Param
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import qualified Data.Array.Knead.Code as Code
 import Data.Array.Knead.Expression (Exp, unExp, )
@@ -42,9 +42,9 @@
 import Control.Monad.HT ((<=<), )
 import Control.Applicative (liftA2, )
 import Data.Tuple.HT (mapFst, )
-import Data.Word (Word32, )
 
 
+{-# INLINE feed #-}
 feed ::
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     MultiValueMemory.C a) =>
@@ -93,13 +93,13 @@
 
 
 foreign import ccall safe "dynamic" callShaper ::
-   Importer (Ptr param -> Ptr shape -> IO Word32)
+   Importer (Ptr param -> Ptr shape -> IO Shape.Size)
 
 
 renderShape ::
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Storable a, MultiValueMemory.C a) =>
-   Sym.Array p sh a -> IO (p -> IO (sh, Word32))
+   Sym.Array p sh a -> IO (p -> IO (sh, Shape.Size))
 renderShape (Sym.Array arr create delete) = do
    fsh <-
       compile "renderShape" $
@@ -129,7 +129,7 @@
 render = PhysHull.render . Sym.arrayHull
 
 
-mapAccumL ::
+mapAccumLSimple ::
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Shape.C n, Storable n, MultiValueMemory.C n,
     MultiValue.C acc,
@@ -139,9 +139,9 @@
    Sym.Array p sh acc ->
    Sym.Array p (sh, n) a ->
    IO (p -> IO (Phys.Array (sh,n) b))
-mapAccumL f arrInit arrMap =
-   PhysHull.mapAccumL $
-      liftA2 (PhysHull.MapAccumL f)
+mapAccumLSimple f arrInit arrMap =
+   PhysHull.mapAccumLSimple $
+      liftA2 (PhysHull.MapAccumLSimple f)
          (Sym.arrayHull arrInit)
          (Sym.arrayHull arrMap)
 
diff --git a/src/Data/Array/Knead/Parameterized/PhysicalHull.hs b/src/Data/Array/Knead/Parameterized/PhysicalHull.hs
--- a/src/Data/Array/Knead/Parameterized/PhysicalHull.hs
+++ b/src/Data/Array/Knead/Parameterized/PhysicalHull.hs
@@ -2,21 +2,32 @@
 {-# LANGUAGE ForeignFunctionInterface #-}
 module Data.Array.Knead.Parameterized.PhysicalHull (
    render,
+   MapFilter(..),
+   mapFilter,
+   FilterOuter(..),
+   filterOuter,
    Scatter(..),
    scatter,
    ScatterMaybe(..),
    scatterMaybe,
+   MapAccumLSimple(..),
+   mapAccumLSimple,
+   MapAccumLSequence(..),
+   mapAccumLSequence,
    MapAccumL(..),
    mapAccumL,
    FoldOuterL(..),
    foldOuterL,
+   AddDimension(..),
+   addDimension,
    ) where
 
 import qualified Data.Array.Knead.Parameterized.Private as Sym
 import qualified Data.Array.Knead.Simple.PhysicalPrivate as Priv
 import qualified Data.Array.Knead.Simple.Physical as Phys
 import qualified Data.Array.Knead.Simple.Private as Core
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
+import qualified Data.Array.Knead.Expression as Expr
 import qualified Data.Array.Knead.Code as Code
 import Data.Array.Knead.Expression (Exp, unExp, )
 import Data.Array.Knead.Code (compile, )
@@ -29,26 +40,41 @@
 import qualified LLVM.Core as LLVM
 
 import Foreign.Marshal.Utils (with, )
+import Foreign.Marshal.Array (allocaArray, )
 import Foreign.Marshal.Alloc (alloca, )
-import Foreign.Storable (Storable, peek, )
-import Foreign.ForeignPtr (withForeignPtr, mallocForeignPtrArray, )
+import Foreign.Storable (Storable, peek, peekElemOff, )
+import Foreign.ForeignPtr (ForeignPtr, withForeignPtr, mallocForeignPtrArray, )
 import Foreign.Ptr (FunPtr, Ptr, )
 
 import Control.Exception (bracket, )
 import Control.Monad.HT (void, )
 import Control.Applicative (liftA2, )
-import Data.Word (Word32, )
 
 
+mallocArray :: (Storable a) => Shape.Size -> IO (ForeignPtr a)
+mallocArray = mallocForeignPtrArray . fromIntegral
+
+withForeignMemPtr ::
+   ForeignPtr a -> (Ptr (MultiValueMemory.Struct a) -> IO b) -> IO b
+withForeignMemPtr fptr act =
+   withForeignPtr fptr $ act . MultiValueMemory.castStructPtr
+
+
 type Importer f = FunPtr f -> f
 
 
 foreign import ccall safe "dynamic" callShaper ::
-   Importer (Ptr param -> Ptr shape -> IO Word32)
+   Importer (Ptr param -> Ptr shape -> IO Shape.Size)
 
 foreign import ccall safe "dynamic" callFill ::
    Importer (Ptr param -> Ptr shape -> Ptr am -> IO ())
 
+
+{-
+Attention:
+The 'fill' function may alter the shape.
+An example is 'mapFilter'.
+-}
 materialize ::
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Storable a, MultiValueMemory.C a) =>
@@ -81,14 +107,121 @@
       with param $ \paramPtr -> do
          let paramMVPtr = MultiValueMemory.castStructPtr paramPtr
          let shapeMVPtr = MultiValueMemory.castStructPtr shptr
-         n <- fsh paramMVPtr shapeMVPtr
-         fptr <- mallocForeignPtrArray (fromIntegral n)
-         withForeignPtr fptr $
-            farr paramMVPtr shapeMVPtr . MultiValueMemory.castStructPtr
+         fptr <- mallocArray =<< fsh paramMVPtr shapeMVPtr
+         withForeignMemPtr fptr $ farr paramMVPtr shapeMVPtr
          sh <- peek shptr
          return (Phys.Array sh fptr)
 
 
+foreign import ccall safe "dynamic" callFillExpArray ::
+   Importer (Ptr param -> Ptr final -> Ptr shape -> Ptr am -> IO ())
+
+
+materializeExpArray ::
+   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
+    Storable a, MultiValueMemory.C a,
+    Storable b, MultiValueMemory.C b) =>
+   String ->
+   (core -> Exp sh) ->
+   (core ->
+    LLVM.Value (Ptr (MultiValueMemory.Struct b)) ->
+    LLVM.Value (Ptr (MultiValueMemory.Struct sh)) ->
+    LLVM.Value (Ptr (MultiValueMemory.Struct a)) ->
+    LLVM.CodeGenFunction () ()) ->
+   Sym.Hull p core -> IO (p -> IO (b, Phys.Array sh a))
+materializeExpArray name shape fill (Sym.Hull core create delete) = do
+   (fsh, farr) <-
+      compile name $
+      liftA2 (,)
+         (Code.createFunction callShaper "shape" $
+          \paramPtr resultPtr -> do
+            param <- Memory.load paramPtr
+            sh <- unExp $ shape $ core param
+            MultiValueMemory.store sh resultPtr
+            Shape.sizeCode sh >>= LLVM.ret)
+         (Code.createFunction callFillExpArray "fill" $
+          \paramPtr finalPtr shapePtr bufferPtr -> do
+            param <- Memory.load paramPtr
+            fill (core param) finalPtr shapePtr bufferPtr
+            LLVM.ret ())
+
+   return $ \p ->
+      bracket (create p) (delete . fst) $ \(_ctx, param) ->
+      alloca $ \shptr ->
+      alloca $ \finalPtr ->
+      with param $ \paramPtr -> do
+         let paramMVPtr = MultiValueMemory.castStructPtr paramPtr
+         let finalMVPtr = MultiValueMemory.castStructPtr finalPtr
+         let shapeMVPtr = MultiValueMemory.castStructPtr shptr
+         fptr <- mallocArray =<< fsh paramMVPtr shapeMVPtr
+         withForeignMemPtr fptr $ farr paramMVPtr finalMVPtr shapeMVPtr
+         sh <- peek shptr
+         final <- peek finalPtr
+         return (final, Phys.Array sh fptr)
+
+
+foreign import ccall safe "dynamic" callShaper2 ::
+   Importer (Ptr param -> Ptr shapeA -> Ptr shapeB -> Ptr Shape.Size -> IO ())
+
+foreign import ccall safe "dynamic" callFill2 ::
+   Importer (Ptr param -> Ptr shapeA -> Ptr am -> Ptr shapeB -> Ptr bm -> IO ())
+
+
+materialize2 ::
+   (Shape.C sha, Storable sha, MultiValueMemory.C sha,
+    Shape.C shb, Storable shb, MultiValueMemory.C shb,
+    Storable a, MultiValueMemory.C a,
+    Storable b, MultiValueMemory.C b) =>
+   String ->
+   (core -> Exp (sha,shb)) ->
+   (core ->
+    (LLVM.Value (Ptr (MultiValueMemory.Struct sha)),
+     LLVM.Value (Ptr (MultiValueMemory.Struct a))) ->
+    (LLVM.Value (Ptr (MultiValueMemory.Struct shb)),
+     LLVM.Value (Ptr (MultiValueMemory.Struct b))) ->
+    LLVM.CodeGenFunction () ()) ->
+   Sym.Hull p core -> IO (p -> IO (Phys.Array sha a, Phys.Array shb b))
+materialize2 name shape fill (Sym.Hull core create delete) = do
+   (fsh, farr) <-
+      compile name $
+      liftA2 (,)
+         (Code.createFunction callShaper2 "shape" $
+          \paramPtr shapeAPtr shapeBPtr sizesPtr -> do
+            param <- Memory.load paramPtr
+            (sha,shb) <- fmap MultiValue.unzip $ unExp $ shape $ core param
+            MultiValueMemory.store sha shapeAPtr
+            MultiValueMemory.store shb shapeBPtr
+            sizeAPtr <- LLVM.bitcast sizesPtr
+            flip LLVM.store sizeAPtr =<< Shape.sizeCode sha
+            sizeBPtr <- A.advanceArrayElementPtr sizeAPtr
+            flip LLVM.store sizeBPtr =<< Shape.sizeCode shb
+            LLVM.ret ())
+         (Code.createFunction callFill2 "fill" $
+          \paramPtr shapeAPtr bufferAPtr shapeBPtr bufferBPtr -> do
+            param <- Memory.load paramPtr
+            fill (core param) (shapeAPtr, bufferAPtr) (shapeBPtr, bufferBPtr)
+            LLVM.ret ())
+
+   return $ \p ->
+      bracket (create p) (delete . fst) $ \(_ctx, param) ->
+      alloca $ \shaPtr ->
+      alloca $ \shbPtr ->
+      allocaArray 2 $ \sizesPtr ->
+      with param $ \paramPtr -> do
+         let paramMVPtr = MultiValueMemory.castStructPtr paramPtr
+         let shapeAMVPtr = MultiValueMemory.castStructPtr shaPtr
+         let shapeBMVPtr = MultiValueMemory.castStructPtr shbPtr
+         fsh paramMVPtr shapeAMVPtr shapeBMVPtr sizesPtr
+         afptr <- mallocArray =<< peekElemOff sizesPtr 0
+         bfptr <- mallocArray =<< peekElemOff sizesPtr 1
+         withForeignMemPtr afptr $ \aptr ->
+            withForeignMemPtr bfptr $ \bptr ->
+            farr paramMVPtr shapeAMVPtr aptr shapeBMVPtr bptr
+         sha <- peek shaPtr
+         shb <- peek shbPtr
+         return (Phys.Array sha afptr, Phys.Array shb bfptr)
+
+
 render ::
    (Shape.C sh, Shape.Index sh ~ ix,
     Storable sh, MultiValueMemory.C sh,
@@ -145,40 +278,151 @@
          Priv.scatterMaybe accum arrInit arrMap)
 
 
-data MapAccumL sh n acc a b =
+data MapAccumLSimple sh n acc a b =
+   MapAccumLSimple {
+      mapAccumLSimpleAccum :: Exp acc -> Exp a -> Exp (acc,b),
+      mapAccumLSimpleInit :: Core.Array sh acc,
+      mapAccumLSimpleArray :: Core.Array (sh, n) a
+   }
+
+mapAccumLSimple ::
+   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
+    Shape.C n, Storable n, MultiValueMemory.C n,
+    MultiValue.C acc,
+    Storable a, MultiValueMemory.C a,
+    Storable b, MultiValueMemory.C b) =>
+   Sym.Hull p (MapAccumLSimple sh n acc a b) -> IO (p -> IO (Phys.Array (sh,n) b))
+mapAccumLSimple =
+   materialize "mapAccumLSimple"
+      (Core.shape . mapAccumLSimpleArray)
+      (\(MapAccumLSimple f arrInit arrData) ->
+         Priv.mapAccumLSimple f arrInit arrData)
+
+
+data MapAccumLSequence n acc final a b =
+   MapAccumLSequence {
+      mapAccumLSequenceAccum :: Exp acc -> Exp a -> Exp (acc,b),
+      mapAccumLSequenceFinal :: Exp acc -> Exp final,
+      mapAccumLSequenceInit :: Exp acc,
+      mapAccumLSequenceArray :: Core.Array n a
+   }
+
+-- FIXME: check correct size of array of initial values
+mapAccumLSequence ::
+   (Shape.C n, Storable n, MultiValueMemory.C n,
+    MultiValue.C acc,
+    Storable final, MultiValueMemory.C final,
+    Storable a, MultiValueMemory.C a,
+    Storable b, MultiValueMemory.C b) =>
+   Sym.Hull p (MapAccumLSequence n acc final a b) ->
+   IO (p -> IO (final, Phys.Array n b))
+mapAccumLSequence =
+   materializeExpArray "mapAccumLSequence"
+      (Core.shape . mapAccumLSequenceArray)
+      (\(MapAccumLSequence f final expInit arr) ->
+         Priv.mapAccumLSequence f final expInit arr)
+
+
+data MapAccumL sh n acc final a b =
    MapAccumL {
       mapAccumLAccum :: Exp acc -> Exp a -> Exp (acc,b),
+      mapAccumLFinal :: Exp acc -> Exp final,
       mapAccumLInit :: Core.Array sh acc,
-      mapAccumLMap :: Core.Array (sh, n) a
+      mapAccumLArray :: Core.Array (sh, n) a
    }
 
+-- FIXME: check correct size of array of initial values
 mapAccumL ::
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Shape.C n, Storable n, MultiValueMemory.C n,
     MultiValue.C acc,
+    Storable final, MultiValueMemory.C final,
     Storable a, MultiValueMemory.C a,
     Storable b, MultiValueMemory.C b) =>
-   Sym.Hull p (MapAccumL sh n acc a b) -> IO (p -> IO (Phys.Array (sh,n) b))
+   Sym.Hull p (MapAccumL sh n acc final a b) ->
+   IO (p -> IO (Phys.Array sh final, Phys.Array (sh,n) b))
 mapAccumL =
-   materialize "mapAccumL"
-      (Core.shape . mapAccumLMap)
-      (\(MapAccumL f arrInit arrData) -> Priv.mapAccumL f arrInit arrData)
+   materialize2 "mapAccumL"
+      (\core ->
+         Expr.zip
+            (Core.shape $ mapAccumLInit core)
+            (Core.shape $ mapAccumLArray core))
+      (\(MapAccumL f final arrInit arrData) ->
+         Priv.mapAccumL f final arrInit arrData)
 
 
 data FoldOuterL n sh a b =
    FoldOuterL {
       foldOuterLAccum :: Exp a -> Exp b -> Exp a,
       foldOuterLInit :: Core.Array sh a,
-      foldOuterLMap :: Core.Array (n,sh) b
+      foldOuterLArray :: Core.Array (n,sh) b
    }
 
 -- FIXME: check correct size of array of initial values
 foldOuterL ::
-   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
-    Shape.C n, Storable n, MultiValueMemory.C n,
+   (Shape.C n, Storable n, MultiValueMemory.C n,
+    Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Storable a, MultiValueMemory.C a) =>
    Sym.Hull p (FoldOuterL n sh a b) -> IO (p -> IO (Phys.Array sh a))
 foldOuterL =
    materialize "foldOuterL"
       (Core.shape . foldOuterLInit)
       (\(FoldOuterL f arrInit arrData) -> Priv.foldOuterL f arrInit arrData)
+
+
+data MapFilter n a b =
+   MapFilter {
+      mapFilterMap :: Exp a -> Exp b,
+      mapFilterPredicate :: Exp a -> Exp Bool,
+      mapFilterArray :: Core.Array n a
+   }
+
+mapFilter ::
+   (Shape.Sequence n,
+    Storable n, MultiValueMemory.C n,
+    Storable b, MultiValueMemory.C b) =>
+   Sym.Hull p (MapFilter n a b) -> IO (p -> IO (Phys.Array n b))
+mapFilter =
+   materialize "mapFilter"
+      (Core.shape . mapFilterArray)
+      (\(MapFilter f p arr) shapePtr bufferPtr ->
+         flip MultiValueMemory.store shapePtr
+            =<< Priv.mapFilter f p arr shapePtr bufferPtr)
+
+
+data FilterOuter n sh a =
+   FilterOuter {
+      filterOuterPredicate :: Core.Array n Bool,
+      filterOuterArray :: Core.Array (n,sh) a
+   }
+
+-- FIXME: check correct size of row selection array
+filterOuter ::
+   (Shape.Sequence n, Storable n, MultiValueMemory.C n,
+    Shape.C sh, Storable sh, MultiValueMemory.C sh,
+    Storable a, MultiValueMemory.C a) =>
+   Sym.Hull p (FilterOuter n sh a) -> IO (p -> IO (Phys.Array (n,sh) a))
+filterOuter =
+   materialize "filterOuter"
+      (Core.shape . filterOuterArray)
+      (\(FilterOuter p arr) shapePtr bufferPtr ->
+         flip MultiValueMemory.store shapePtr
+            =<< Priv.filterOuter p arr shapePtr bufferPtr)
+
+
+data AddDimension sh n a b =
+   AddDimension {
+      addDimensionSize :: Exp n,
+      addDimensionSelect :: Exp (Shape.Index n) -> Exp a -> Exp b,
+      addDimensionArray :: Core.Array sh a
+   }
+
+addDimension ::
+   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
+    Shape.C n, Storable n, MultiValueMemory.C n,
+    Storable b, MultiValueMemory.C b) =>
+   Sym.Hull p (AddDimension sh n a b) -> IO (p -> IO (Phys.Array (sh,n) b))
+addDimension =
+   materialize "addDimension"
+      (\r -> Expr.zip (Core.shape (addDimensionArray r)) (addDimensionSize r))
+      (\(AddDimension n select arr) -> Priv.addDimension n select arr)
diff --git a/src/Data/Array/Knead/Parameterized/Private.hs b/src/Data/Array/Knead/Parameterized/Private.hs
--- a/src/Data/Array/Knead/Parameterized/Private.hs
+++ b/src/Data/Array/Knead/Parameterized/Private.hs
@@ -5,7 +5,7 @@
 import qualified Data.Array.Knead.Simple.Symbolic as Core
 
 import qualified Data.Array.Knead.Parameter as Param
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp, )
 
@@ -113,7 +113,7 @@
     MultiValueMemory.C c, Storable c, MultiValue.C a) =>
    (Exp c -> Exp a -> Exp a -> Exp a) ->
    Param.T p c -> Array p sh a -> Array p z a
-fold1All = lift Core.fold1All
+fold1All = lift (\p -> Core.fill Shape.scalar . Core.fold1All p)
 
 lift ::
    (Shape.C sh0, Shape.C sh1,
@@ -200,6 +200,7 @@
 deletePlain () = return ()
 
 
+{-# INLINE combineCreate #-}
 combineCreate ::
    Monad m =>
    (p -> m (ctxA, paramA)) -> (p -> m (ctxB, paramB)) ->
@@ -209,6 +210,7 @@
    (ctxB, paramB) <- createB p
    return ((ctxA, ctxB), (paramA, paramB))
 
+{-# INLINE combineDelete #-}
 combineDelete ::
    Monad m =>
    (ctxA -> m ()) -> (ctxB -> m ()) -> (ctxA, ctxB) -> m ()
diff --git a/src/Data/Array/Knead/Parameterized/Render.hs b/src/Data/Array/Knead/Parameterized/Render.hs
--- a/src/Data/Array/Knead/Parameterized/Render.hs
+++ b/src/Data/Array/Knead/Parameterized/Render.hs
@@ -4,10 +4,15 @@
 -}
 module Data.Array.Knead.Parameterized.Render (
    run,
+   MapFilter(..),
+   FilterOuter(..),
    Scatter(..),
    ScatterMaybe(..),
+   MapAccumLSimple(..),
+   MapAccumLSequence(..),
    MapAccumL(..),
    FoldOuterL(..),
+   AddDimension(..),
    ) where
 
 import qualified Data.Array.Knead.Parameterized.PhysicalHull as PhysHullP
@@ -16,9 +21,11 @@
 import qualified Data.Array.Knead.Simple.Physical as Phys
 import qualified Data.Array.Knead.Simple.Private as Core
 import qualified Data.Array.Knead.Parameter as Param
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import Data.Array.Knead.Parameterized.PhysicalHull
-         (Scatter, ScatterMaybe, MapAccumL, FoldOuterL)
+         (MapFilter, FilterOuter,
+          MapAccumLSimple, MapAccumLSequence, MapAccumL, FoldOuterL,
+          Scatter, ScatterMaybe, AddDimension)
 import Data.Array.Knead.Expression (Exp, )
 
 import qualified LLVM.Extra.Multi.Value.Memory as MultiValueMemory
@@ -44,6 +51,22 @@
    build = PhysHullP.render
 
 instance
+   (Shape.Sequence n,
+    Storable n, MultiValueMemory.C n,
+    Storable b, MultiValueMemory.C b) =>
+      C (MapFilter n a b) where
+   type Plain (MapFilter n a b) = IO (Phys.Array n b)
+   build = PhysHullP.mapFilter
+
+instance
+   (Shape.Sequence n, Storable n, MultiValueMemory.C n,
+    Shape.C sh, Storable sh, MultiValueMemory.C sh,
+    Storable a, MultiValueMemory.C a) =>
+      C (FilterOuter n sh a) where
+   type Plain (FilterOuter n sh a) = IO (Phys.Array (n,sh) a)
+   build = PhysHullP.filterOuter
+
+instance
    (MultiValueMemory.C sh0, Storable sh0, Shape.C sh0,
     MultiValueMemory.C sh1, Storable sh1, Shape.C sh1,
     MultiValueMemory.C a, Storable a) =>
@@ -65,18 +88,48 @@
     MultiValue.C acc,
     Storable a, MultiValueMemory.C a,
     Storable b, MultiValueMemory.C b) =>
-      C (MapAccumL sh n acc a b) where
-   type Plain (MapAccumL sh n acc a b) = IO (Phys.Array (sh,n) b)
-   build = PhysHullP.mapAccumL
+      C (MapAccumLSimple sh n acc a b) where
+   type Plain (MapAccumLSimple sh n acc a b) = IO (Phys.Array (sh,n) b)
+   build = PhysHullP.mapAccumLSimple
 
 instance
+   (Shape.C n, Storable n, MultiValueMemory.C n,
+    MultiValue.C acc,
+    Storable final, MultiValueMemory.C final,
+    Storable a, MultiValueMemory.C a,
+    Storable b, MultiValueMemory.C b) =>
+      C (MapAccumLSequence n acc final a b) where
+   type Plain (MapAccumLSequence n acc final a b) = IO (final, Phys.Array n b)
+   build = PhysHullP.mapAccumLSequence
+
+instance
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Shape.C n, Storable n, MultiValueMemory.C n,
+    MultiValue.C acc,
+    Storable final, MultiValueMemory.C final,
     Storable a, MultiValueMemory.C a,
     Storable b, MultiValueMemory.C b) =>
+      C (MapAccumL sh n acc final a b) where
+   type Plain (MapAccumL sh n acc final a b) =
+            IO (Phys.Array sh final, Phys.Array (sh,n) b)
+   build = PhysHullP.mapAccumL
+
+instance
+   (Shape.C n, Storable n, MultiValueMemory.C n,
+    Shape.C sh, Storable sh, MultiValueMemory.C sh,
+    Storable a, MultiValueMemory.C a,
+    Storable b, MultiValueMemory.C b) =>
       C (FoldOuterL n sh a b) where
    type Plain (FoldOuterL n sh a b) = IO (Phys.Array sh a)
    build = PhysHullP.foldOuterL
+
+instance
+   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
+    Shape.C n, Storable n, MultiValueMemory.C n,
+    Storable b, MultiValueMemory.C b) =>
+      C (AddDimension sh n a b) where
+   type Plain (AddDimension sh n a b) = IO (Phys.Array (sh,n) b)
+   build = PhysHullP.addDimension
 
 
 singleton :: Exp a -> Core.Array () a
diff --git a/src/Data/Array/Knead/Parameterized/Slice.hs b/src/Data/Array/Knead/Parameterized/Slice.hs
--- a/src/Data/Array/Knead/Parameterized/Slice.hs
+++ b/src/Data/Array/Knead/Parameterized/Slice.hs
@@ -17,12 +17,12 @@
 import qualified Data.Array.Knead.Simple.Slice as Slice
 import qualified Data.Array.Knead.Simple.Private as Core
 
-import qualified Data.Array.Knead.Index.Linear as Linear
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Cubic as Linear
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Parameter as Param
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp, )
-import Data.Array.Knead.Index.Linear ((:.), )
+import Data.Array.Knead.Shape.Cubic ((:.), )
 
 import qualified LLVM.Extra.Multi.Value.Memory as MultiValueMemory
 import qualified LLVM.Extra.Multi.Value as MultiValue
diff --git a/src/Data/Array/Knead/Parameterized/Symbolic.hs b/src/Data/Array/Knead/Parameterized/Symbolic.hs
--- a/src/Data/Array/Knead/Parameterized/Symbolic.hs
+++ b/src/Data/Array/Knead/Parameterized/Symbolic.hs
@@ -22,7 +22,7 @@
 import qualified Data.Array.Knead.Simple.Symbolic as Core
 import Data.Array.Knead.Parameterized.Private (Array, gather, )
 
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Parameter as Param
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp, )
diff --git a/src/Data/Array/Knead/Shape/Cubic.hs b/src/Data/Array/Knead/Shape/Cubic.hs
new file mode 100644
--- /dev/null
+++ b/src/Data/Array/Knead/Shape/Cubic.hs
@@ -0,0 +1,581 @@
+{-# LANGUAGE Rank2Types #-}
+{-# LANGUAGE ExistentialQuantification #-}
+{-# LANGUAGE EmptyDataDecls #-}
+{-# LANGUAGE TypeFamilies #-}
+{-# LANGUAGE TypeOperators #-}
+{-# LANGUAGE FlexibleContexts #-}
+{-# LANGUAGE FlexibleInstances #-}
+{-# LANGUAGE UndecidableInstances #-}
+{-# LANGUAGE MultiParamTypeClasses #-}
+module Data.Array.Knead.Shape.Cubic (
+   C(switch),
+   switchInt,
+   intersect,
+   value,
+   constant,
+   paramWith,
+   tunnel,
+   flattenIndex,
+   peek,
+   poke,
+   computeSize,
+
+   Struct,
+   T(..),
+   Z(Z), z,
+   (:.)((:.)),
+   Shape, shape,
+   Index, index,
+   cons, (#:.),
+   head,
+   tail,
+   switchR,
+   loadMultiValue,
+   storeMultiValue,
+   ) where
+
+import qualified Data.Array.Knead.Shape.Nested as Shape
+import qualified Data.Array.Knead.Shape.Cubic.Int as Index
+
+import qualified Data.Array.Knead.Parameter as Param
+import qualified Data.Array.Knead.Expression as Expr
+import Data.Array.Knead.Expression (Exp, )
+
+import qualified LLVM.Extra.Multi.Value.Memory as MultiValueMemory
+import qualified LLVM.Extra.Multi.Value as MultiValue
+import qualified LLVM.Extra.Multi.Iterator as IterMV
+import qualified LLVM.Extra.Iterator as Iter
+import qualified LLVM.Extra.Arithmetic as A
+import qualified LLVM.Extra.Control as C
+import LLVM.Extra.Multi.Value (Atom, )
+
+import qualified LLVM.Util.Loop as Loop
+import qualified LLVM.Core as LLVM
+
+import qualified Foreign.Storable as St
+import Foreign.Storable.FixedArray (sizeOfArray, )
+import Foreign.Marshal.Array (advancePtr, )
+import Foreign.Ptr (Ptr, castPtr, )
+
+import Control.Monad (liftM2, )
+
+import Prelude hiding (min, head, tail, )
+
+
+class C ix where
+   switch ::
+      f Z ->
+      (forall ix0 i. (C ix0, Index.Single i) => f (ix0 :. i)) ->
+      f ix
+
+instance C Z where
+   switch x _ = x
+
+instance (C ix0, Index.Single i) => C (ix0 :. i) where
+   switch _ x = x
+
+
+newtype SwitchInt f ix i = SwitchInt {runSwitchInt :: f (ix :. i)}
+
+switchInt ::
+   (C ix) =>
+   f Z ->
+   (forall ix0. (C ix0) => f (ix0 :. Index.Int)) ->
+   f ix
+switchInt z0 cons0 =
+   switch z0
+      (runSwitchInt $ Index.switchSingle (SwitchInt cons0))
+
+
+newtype Op2 tag sh = Op2 {runOp2 :: Exp (T tag sh) -> Exp (T tag sh) -> Exp (T tag sh)}
+
+intersect :: C sh => Exp (Shape sh) -> Exp (Shape sh) -> Exp (Shape sh)
+intersect =
+   runOp2 $
+   switchInt
+      (Op2 $ \z0 _ -> z0)
+      (Op2 $
+       switchR $ \is i ->
+       switchR $ \js j ->
+          intersect is js #:. Expr.min i j)
+
+
+_value :: (C sh, MultiValue.C sh) => sh -> Exp sh
+_value = Expr.lift0 . MultiValue.cons
+
+
+newtype MakeValue val tag sh = MakeValue {runMakeValue :: T tag sh -> val (T tag sh)}
+
+value :: (C sh, Expr.Value val) => T tag sh -> val (T tag sh)
+value =
+   runMakeValue $
+   switchInt
+      (MakeValue $ \(Cons Z) -> z)
+      (MakeValue $ \(Cons (t:.h)) ->
+         value (Cons t) #:. Expr.lift0 (MultiValue.cons h))
+
+paramWith ::
+   (C sh, Expr.Value val) =>
+   Param.T p (T tag sh) ->
+   (forall parameters.
+    (St.Storable parameters,
+     MultiValueMemory.C parameters) =>
+    (p -> parameters) ->
+    (MultiValue.T parameters -> val (T tag sh)) ->
+    a) ->
+   a
+paramWith p f =
+   case tunnel p of
+      Param.Tunnel get val -> f get (Expr.lift0 . val)
+
+tunnel :: (C sh) => Param.T p (T tag sh) -> Param.Tunnel p (T tag sh)
+tunnel p =
+   case structFieldsPropF p of
+      StructFieldsProp -> Param.tunnel value p
+
+
+data StructFieldsProp sh = LLVM.StructFields (Struct sh) => StructFieldsProp
+
+_structFieldsProp :: (C sh) => f sh -> StructFieldsProp sh
+_structFieldsProp _p = structFieldsRec
+
+structFieldsPropF :: (C sh) => f (g sh) -> StructFieldsProp sh
+structFieldsPropF _p = structFieldsRec
+
+withStructFieldsPropFF ::
+   (C sh) => (StructFieldsProp sh -> f (g (h sh))) -> f (g (h sh))
+withStructFieldsPropFF f = f structFieldsRec
+
+structFieldsRec :: (C sh) => StructFieldsProp sh
+structFieldsRec =
+   switchInt
+      StructFieldsProp
+      (succStructFieldsProp structFieldsRec)
+
+succStructFieldsProp ::
+   StructFieldsProp sh -> StructFieldsProp (sh:.Index.Int)
+succStructFieldsProp StructFieldsProp = StructFieldsProp
+
+
+data Z = Z
+   deriving (Eq, Ord, Read, Show)
+
+
+infixl 3 :., #:.
+
+data tail :. head = !tail :. !head
+   deriving (Eq, Ord, Read, Show)
+
+
+newtype T tag sh = Cons {decons :: sh}
+
+data ShapeTag
+data IndexTag
+
+type Shape = T ShapeTag
+type Index = T IndexTag
+
+shape :: sh -> Shape sh
+shape = Cons
+
+index :: ix -> Index ix
+index = Cons
+
+
+(#:.) :: (Expr.Value val) => val (T tag sh) -> val i -> val (T tag (sh:.i))
+(#:.) = cons
+
+cons :: (Expr.Value val) => val (T tag sh) -> val i -> val (T tag (sh:.i))
+cons =
+   Expr.lift2 $
+      \(MultiValue.Cons t) (MultiValue.Cons h) ->
+         MultiValue.Cons (t,h)
+
+z :: (Expr.Value val) => val (T tag Z)
+z = Expr.lift0 $ MultiValue.Cons ()
+
+head :: (Expr.Value val) => val (T tag (sh:.i)) -> val i
+head = Expr.lift1 $ \(MultiValue.Cons (_t,h)) -> MultiValue.Cons h
+
+tail :: (Expr.Value val) => val (T tag (sh:.i)) -> val (T tag sh)
+tail = Expr.lift1 $ \(MultiValue.Cons (t,_h)) -> MultiValue.Cons t
+
+switchR ::
+   Expr.Value val =>
+   (val (T tag sh) -> val i -> a) -> val (T tag (sh :. i)) -> a
+switchR f ix = f (tail ix) (head ix)
+
+
+instance (tag ~ ShapeTag, sh ~ Z) => Shape.Scalar (T tag sh) where
+   scalar = Expr.lift0 $ MultiValue.Cons ()
+   zeroIndex _ = Expr.lift0 $ MultiValue.Cons ()
+
+
+type family PatternTuple pattern
+type family Decomposed (f :: * -> *) tag pattern
+
+type instance PatternTuple (sh:.s) =
+   PatternTuple sh :. MultiValue.PatternTuple s
+
+type instance Decomposed f tag (sh:.s) =
+   Decomposed f tag sh :. MultiValue.Decomposed f s
+
+type instance PatternTuple (Atom sh) = sh
+
+type instance Decomposed f tag (Atom sh) = f (T tag sh)
+
+
+class
+   (Expr.Composed (Decomposed Exp tag pattern) ~ T tag (PatternTuple pattern)) =>
+      Decompose tag pattern where
+   decompose ::
+      T tag pattern -> Exp (T tag (PatternTuple pattern)) ->
+      Decomposed Exp tag pattern
+
+instance Decompose tag (Atom sh) where
+   decompose (Cons _atom) x = x
+
+instance (Decompose tag sh, Expr.Decompose s) => Decompose tag (sh :. s) where
+   decompose (Cons (psh:.ps)) x =
+      decompose (Cons psh) (tail x) :. Expr.decompose ps (head x)
+
+
+type instance MultiValue.PatternTuple (T tag sh) = T tag (PatternTuple sh)
+
+type instance MultiValue.Decomposed f (T tag sh) = Decomposed f tag sh
+
+
+type family Unwrap sh
+type instance Unwrap (T tag sh) = sh
+
+type family Tag sh
+type instance Tag (T tag sh) = tag
+
+instance
+   (Expr.Compose sh,
+    Expr.Composed sh ~ T (Tag (Expr.Composed sh)) (Unwrap (Expr.Composed sh)),
+    Expr.Compose s) =>
+      Expr.Compose (sh :. s) where
+   type Composed (sh :. s) =
+           T (Tag (Expr.Composed sh))
+             (Unwrap (Expr.Composed sh) :. Expr.Composed s)
+   compose (sh :. s) = cons (Expr.compose sh) (Expr.compose s)
+
+instance (Decompose tag sh) => Expr.Decompose (T tag sh) where
+   decompose = decompose
+
+
+
+instance (C sh) => St.Storable (T tag sh) where
+   sizeOf (Cons sh) = sizeOfArray (rank sh) (0::Shape.Size)
+   alignment (Cons _sh) = St.alignment (0::Shape.Size)
+   poke ptr = poke (castPtr ptr) . decons
+   peek = fmap Cons . peek . castPtr
+
+
+type family Repr (f :: * -> *) sh
+type instance Repr f Z = ()
+type instance Repr f (tail :. head) = (Repr f tail, MultiValue.Repr f head)
+
+instance (C sh) => MultiValue.C (T tag sh) where
+   type Repr f (T tag sh) = Repr f sh
+   cons = value
+   undef = constant $ MultiValue.undef
+   zero = constant $ MultiValue.zero
+   addPhis = addPhis
+   phis = phis
+
+instance (tag ~ ShapeTag, C sh) => Shape.C (T tag sh) where
+   type Index (T tag sh) = Index sh
+   size = fromIntegral . size . decons
+   sizeCode = computeSize
+   intersectCode = Expr.unliftM2 intersect
+   flattenIndexRec sh ix =
+      -- a joint implementation would not be more efficient
+      liftM2 (,)
+         (computeSize sh)
+         (flattenIndex sh ix)
+   iterator = iterator
+   loop = loop
+
+
+type family Struct sh
+type instance Struct Z = ()
+type instance Struct (sh :. Index.Int) = (Shape.Size, Struct sh)
+
+instance
+   (C sh, LLVM.StructFields (Struct sh)) =>
+      MultiValueMemory.C (T tag sh) where
+   type Struct (T tag sh) = LLVM.Struct (Struct sh)
+   load = loadMultiValue
+   store = storeMultiValue
+
+loadMultiValue ::
+   (C sh) =>
+   LLVM.Value (Ptr (LLVM.Struct (Struct sh))) ->
+   LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
+loadMultiValue ptr =
+   withStructFieldsPropFF $ \StructFieldsProp ->
+      load =<< castPtrValue ptr
+
+storeMultiValue ::
+   (C sh) =>
+   MultiValue.T (T tag sh) ->
+   LLVM.Value (Ptr (LLVM.Struct (Struct sh))) -> LLVM.CodeGenFunction r ()
+storeMultiValue x ptr =
+   case structFieldsPropF x of
+      StructFieldsProp -> store x =<< castPtrValue ptr
+
+
+newtype FlattenIndex r sh =
+   FlattenIndex {
+      runFlattenIndex ::
+         MultiValue.T (Shape sh) -> MultiValue.T (Index sh) ->
+         LLVM.CodeGenFunction r (LLVM.Value Shape.Size)
+   }
+
+flattenIndex ::
+   (C sh) =>
+   MultiValue.T (Shape sh) -> MultiValue.T (Index sh) ->
+   LLVM.CodeGenFunction r (LLVM.Value Shape.Size)
+flattenIndex =
+   runFlattenIndex $
+   switchInt
+      (FlattenIndex $ \_zerosh _zeroix -> return A.zero)
+      (FlattenIndex $
+         switchR $ \sh (MultiValue.Cons s) ->
+         switchR $ \ix (MultiValue.Cons i) ->
+            A.add i =<< A.mul s =<< flattenIndex sh ix)
+
+
+newtype Rank sh = Rank {runRank :: sh -> Int}
+
+rank :: (C sh) => sh -> Int
+rank =
+   runRank $
+   switch
+      (Rank $ const 0)
+      (Rank $ succ . rank . (\(sh :. _s) -> sh))
+
+
+newtype Peek sh = Peek {runPeek :: Ptr Shape.Size -> IO sh}
+
+peek :: (C sh) => Ptr Shape.Size -> IO sh
+peek =
+   runPeek $
+   switchInt
+      (Peek $ const $ return Z)
+      (Peek $ \ptr -> do
+         h <- St.peek ptr
+         t <- peek $ advancePtr ptr 1
+         return (t :. Index.Int h))
+
+
+newtype Poke sh = Poke {runPoke :: Ptr Shape.Size -> sh -> IO ()}
+
+poke :: (C sh) => Ptr Shape.Size -> sh -> IO ()
+poke =
+   runPoke $
+   switchInt
+      (Poke $ const $ const $ return ())
+      (Poke $ \ptr (sh :. Index.Int i) -> do
+         St.poke ptr i
+         poke (advancePtr ptr 1) sh)
+
+
+castPtrValue ::
+   (LLVM.StructFields sh) =>
+   LLVM.Value (Ptr (LLVM.Struct sh)) ->
+   LLVM.CodeGenFunction r (LLVM.Value (Ptr Shape.Size))
+castPtrValue = LLVM.bitcast
+
+newtype Load r tag sh =
+   Load {
+      runLoad ::
+         LLVM.Value (Ptr Shape.Size) ->
+         LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
+   }
+
+load ::
+   (C sh) =>
+   LLVM.Value (Ptr Shape.Size) ->
+   LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
+load =
+   runLoad $
+   switchInt
+      (Load $ const $ return z)
+      (Load $ \ptr -> do
+         h <- LLVM.load ptr
+         t <- load =<< A.advanceArrayElementPtr ptr
+         return (t #:. MultiValue.Cons h))
+
+
+newtype Store r tag sh =
+   Store {
+      runStore ::
+         MultiValue.T (T tag sh) ->
+         LLVM.Value (Ptr Shape.Size) ->
+         LLVM.CodeGenFunction r ()
+   }
+
+store ::
+   (C sh) =>
+   MultiValue.T (T tag sh) ->
+   LLVM.Value (Ptr Shape.Size) ->
+   LLVM.CodeGenFunction r ()
+store =
+   runStore $
+   switchInt
+      (Store $ \_z _ptr -> return ())
+      (Store $ switchR $ \sh (MultiValue.Cons k) ptr -> do
+         LLVM.store k ptr
+         store sh =<< A.advanceArrayElementPtr ptr)
+
+
+newtype Size sh =
+   Size {
+      runSize :: sh -> Shape.Size
+   }
+
+size :: (C sh) => sh -> Shape.Size
+size =
+   runSize $
+   switchInt
+      (Size $ \_z -> 1)
+      (Size $ \(sh :. Index.Int k) -> k * size sh)
+
+
+newtype ComputeSize r sh =
+   ComputeSize {
+      runComputeSize ::
+         MultiValue.T (Shape sh) ->
+         LLVM.CodeGenFunction r (LLVM.Value Shape.Size)
+   }
+
+computeSize ::
+   (C sh) =>
+   MultiValue.T (Shape sh) ->
+   LLVM.CodeGenFunction r (LLVM.Value Shape.Size)
+computeSize =
+   runComputeSize $
+   switchInt
+      (ComputeSize $ \_z -> return A.one)
+      (ComputeSize $ switchR $ \sh (MultiValue.Cons k) ->
+         A.mul k =<< computeSize sh)
+
+
+newtype
+   Constant val tag sh =
+      Constant {getConstant :: val Index.Int -> val (T tag sh)}
+
+constant :: (C sh, Expr.Value val) => val Index.Int -> val (T tag sh)
+constant =
+   getConstant $
+   switchInt
+      (Constant $ const z)
+      (Constant $ \x -> constant x #:. x)
+
+
+newtype AddPhis r tag sh =
+   AddPhis {
+      runAddPhis ::
+         LLVM.BasicBlock ->
+         MultiValue.T (T tag sh) ->
+         MultiValue.T (T tag sh) ->
+         LLVM.CodeGenFunction r ()
+   }
+
+addPhis ::
+   (C sh) =>
+   LLVM.BasicBlock ->
+   MultiValue.T (T tag sh) ->
+   MultiValue.T (T tag sh) ->
+   LLVM.CodeGenFunction r ()
+addPhis =
+   runAddPhis $
+   switchInt
+      (AddPhis $ \_ _ _ -> return ())
+      (AddPhis $ \bb ->
+       switchR $ \hx tx ->
+       switchR $ \hy ty ->
+          MultiValue.addPhis bb tx ty >>
+          addPhis bb hx hy)
+
+
+newtype Phis r tag sh =
+   Phis {
+      runPhis ::
+         LLVM.BasicBlock ->
+         MultiValue.T (T tag sh) ->
+         LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
+   }
+
+phis ::
+   (C sh) =>
+   LLVM.BasicBlock ->
+   MultiValue.T (T tag sh) ->
+   LLVM.CodeGenFunction r (MultiValue.T (T tag sh))
+phis =
+   runPhis $
+   switchInt
+      (Phis $ \_ -> return)
+      (Phis $ \bb ->
+       switchR $ \h t ->
+          liftM2 (#:.)
+             (phis bb h)
+             (MultiValue.phis bb t))
+
+
+newtype Iterator r sh =
+   Iterator {
+      runIterator ::
+         MultiValue.T (Shape sh) -> Iter.T r (MultiValue.T (Index sh))
+   }
+
+iterator ::
+   (C sh) =>
+   MultiValue.T (Shape sh) -> Iter.T r (MultiValue.T (Index sh))
+iterator =
+   runIterator $
+   switchInt
+      (Iterator $ \ _z -> Iter.empty)
+      (Iterator $ switchR $ \sh n ->
+       fmap (\(ix,i) -> ix#:.i) $
+       Iter.cartesian
+         (iterator sh)
+         (IterMV.takeWhile (MultiValue.cmp LLVM.CmpGT n) $
+          Iter.iterate MultiValue.inc MultiValue.zero))
+
+
+newtype Loop r state sh =
+   Loop {
+      runLoop ::
+         (MultiValue.T (Index sh) ->
+          state ->
+          LLVM.CodeGenFunction r state) ->
+         MultiValue.T (Shape sh) ->
+         state ->
+         LLVM.CodeGenFunction r state
+   }
+
+loop ::
+   (C sh, Loop.Phi state) =>
+   (MultiValue.T (Index sh) ->
+    state ->
+    LLVM.CodeGenFunction r state) ->
+   MultiValue.T (Shape sh) ->
+   state ->
+   LLVM.CodeGenFunction r state
+loop =
+   runLoop $
+   switchInt
+      (Loop $ \code _z -> code z)
+      (Loop $ \code -> switchR $ \sh (MultiValue.Cons n) ->
+         loop
+            (\ix ptrStart ->
+               fmap fst $
+               C.fixedLengthLoop n (ptrStart, A.zero) $ \(ptr, k) ->
+                  liftM2 (,)
+                     (code (ix #:. MultiValue.Cons k) ptr)
+                     (A.inc k))
+            sh)
diff --git a/src/Data/Array/Knead/Shape/Cubic/Int.hs b/src/Data/Array/Knead/Shape/Cubic/Int.hs
new file mode 100644
--- /dev/null
+++ b/src/Data/Array/Knead/Shape/Cubic/Int.hs
@@ -0,0 +1,59 @@
+{-# LANGUAGE TypeFamilies #-}
+module Data.Array.Knead.Shape.Cubic.Int (
+   Single(..),
+   Int(Int), cons, decons,
+   ) where
+
+import qualified Data.Array.Knead.Expression as Expr
+
+import qualified LLVM.Extra.Multi.Value as MultiValue
+import qualified LLVM.Extra.Arithmetic as A
+
+import Data.Word (Word64, )
+
+import Prelude hiding (Int, head, tail, )
+
+
+newtype Int = Int Word64
+
+cons :: (Expr.Value val) => val Word64 -> val Int
+cons = Expr.lift1 $ \(MultiValue.Cons x) -> MultiValue.Cons x
+
+decons :: (Expr.Value val) => val Int -> val Word64
+decons = Expr.lift1 $ \(MultiValue.Cons x) -> MultiValue.Cons x
+
+
+class Single ix where
+   switchSingle :: f Int -> f ix
+
+instance Single Int where
+   switchSingle x = x
+
+
+instance MultiValue.C Int where
+   type Repr f Int = f Word64
+   cons (Int x) = MultiValue.consPrimitive x
+   undef = MultiValue.undefPrimitive
+   zero = MultiValue.zeroPrimitive
+   phis = MultiValue.phisPrimitive
+   addPhis = MultiValue.addPhisPrimitive
+
+instance MultiValue.Additive Int where
+   add = MultiValue.liftM2 A.add
+   sub = MultiValue.liftM2 A.sub
+   neg = MultiValue.liftM A.neg
+
+instance MultiValue.PseudoRing Int where
+   mul = MultiValue.liftM2 A.mul
+
+instance MultiValue.Real Int where
+   min = MultiValue.liftM2 A.min
+   max = MultiValue.liftM2 A.max
+   abs = MultiValue.liftM A.abs
+   signum = MultiValue.liftM A.signum
+
+instance MultiValue.IntegerConstant Int where
+   fromInteger' = cons . A.fromInteger'
+
+instance MultiValue.Comparison Int where
+   cmp mode = MultiValue.liftM2 $ A.cmp mode
diff --git a/src/Data/Array/Knead/Shape/Nested.hs b/src/Data/Array/Knead/Shape/Nested.hs
new file mode 100644
--- /dev/null
+++ b/src/Data/Array/Knead/Shape/Nested.hs
@@ -0,0 +1,539 @@
+{-# LANGUAGE TypeFamilies #-}
+{-# LANGUAGE Rank2Types #-}
+{-# LANGUAGE FlexibleContexts #-}
+module Data.Array.Knead.Shape.Nested (
+   C(..),
+   Size,
+   value,
+   paramWith,
+   load,
+   intersect,
+   flattenIndex,
+
+   Range(..),
+   Shifted(..),
+
+   Scalar(..),
+   Sequence(..),
+   ) where
+
+import qualified Data.Array.Knead.Expression as Expr
+import qualified Data.Array.Knead.Parameter as Param
+import Data.Array.Knead.Expression (Exp, )
+
+import qualified LLVM.Extra.Multi.Value.Memory as MultiMem
+import qualified LLVM.Extra.Multi.Value as MultiValue
+import qualified LLVM.Extra.Multi.Iterator as IterMV
+import qualified LLVM.Extra.Iterator as Iter
+import qualified LLVM.Extra.Arithmetic as A
+import LLVM.Extra.Multi.Value (atom)
+
+import qualified LLVM.Util.Loop as Loop
+import qualified LLVM.Core as LLVM
+
+import qualified Type.Data.Num.Decimal as TypeNum
+
+import Foreign.Storable
+         (Storable, sizeOf, alignment, poke, peek, pokeElemOff, peekElemOff)
+import Foreign.Ptr (Ptr, castPtr)
+
+import Data.Word (Word8, Word16, Word32, Word64)
+import Data.Int (Int8, Int16, Int32, Int64)
+
+import qualified Control.Monad.HT as Monad
+import Control.Applicative ((<$>))
+
+
+type Size = Word64
+
+value :: (C sh, Expr.Value val) => sh -> val sh
+value = Expr.lift0 . MultiValue.cons
+
+paramWith ::
+   (Storable b, MultiMem.C b, Expr.Value val) =>
+   Param.T p b ->
+   (forall parameters.
+    (Storable parameters, MultiMem.C parameters) =>
+    (p -> parameters) ->
+    (MultiValue.T parameters -> val b) ->
+    a) ->
+   a
+paramWith p f =
+   Param.withMulti p (\get val -> f get (Expr.lift0 . val))
+
+load ::
+   (MultiMem.C sh) =>
+   f sh -> LLVM.Value (Ptr (MultiMem.Struct sh)) ->
+   LLVM.CodeGenFunction r (MultiValue.T sh)
+load _ = MultiMem.load
+
+intersect :: (C sh) => Exp sh -> Exp sh -> Exp sh
+intersect = Expr.liftM2 intersectCode
+
+flattenIndex ::
+   (C sh) =>
+   MultiValue.T sh -> MultiValue.T (Index sh) ->
+   LLVM.CodeGenFunction r (LLVM.Value Size)
+flattenIndex sh ix =
+   fmap snd $ flattenIndexRec sh ix
+
+class (MultiValue.C sh, MultiValue.C (Index sh)) => C sh where
+   type Index sh :: *
+   {-
+   It would be better to restrict zipWith to matching shapes
+   and turn shape intersection into a bound check.
+   -}
+   intersectCode ::
+      MultiValue.T sh -> MultiValue.T sh ->
+      LLVM.CodeGenFunction r (MultiValue.T sh)
+   sizeCode ::
+      MultiValue.T sh ->
+      LLVM.CodeGenFunction r (LLVM.Value Size)
+   size :: sh -> Int
+   {- |
+   Result is @(size, flattenedIndex)@.
+   @size@ must equal the result of 'sizeCode'.
+   We use this for sharing intermediate results.
+   -}
+   flattenIndexRec ::
+      MultiValue.T sh -> MultiValue.T (Index sh) ->
+      LLVM.CodeGenFunction r (LLVM.Value Size, LLVM.Value Size)
+   iterator :: (Index sh ~ ix) => MultiValue.T sh -> Iter.T r (MultiValue.T ix)
+   loop ::
+      (Index sh ~ ix, MultiValue.C ix, Loop.Phi state) =>
+      (MultiValue.T ix -> state -> LLVM.CodeGenFunction r state) ->
+      MultiValue.T sh -> state -> LLVM.CodeGenFunction r state
+   loop f sh = Iter.mapState_ f (iterator sh)
+
+
+instance C () where
+   type Index () = ()
+   intersectCode _ _ = return $ MultiValue.cons ()
+   sizeCode _ = return A.one
+   size _ = 1
+   flattenIndexRec _ _ = return (A.one, A.zero)
+   iterator = Iter.singleton
+   loop = id
+
+
+class C sh => Scalar sh where
+   scalar :: (Expr.Value val) => val sh
+   zeroIndex :: (Expr.Value val) => f sh -> val (Index sh)
+
+instance Scalar () where
+   scalar = Expr.lift0 $ MultiValue.Cons ()
+   zeroIndex _ = Expr.lift0 $ MultiValue.Cons ()
+
+
+class
+   (C sh,
+    MultiValue.IntegerConstant (Index sh),
+    MultiValue.Additive (Index sh)) =>
+      Sequence sh where
+   sequenceShapeFromIndex ::
+      MultiValue.T (Index sh) -> LLVM.CodeGenFunction r (MultiValue.T sh)
+
+
+iteratorPrimitive ::
+   (MultiValue.Repr LLVM.Value j ~ LLVM.Value j,
+    Num j, LLVM.IsConst j, LLVM.IsInteger j,
+    LLVM.CmpRet j, LLVM.CmpResult j ~ Bool,
+    MultiValue.Additive j, MultiValue.IntegerConstant j) =>
+   MultiValue.T j -> Iter.T r (MultiValue.T j)
+iteratorPrimitive (MultiValue.Cons n) = iteratorStart n MultiValue.zero
+
+iteratorStart ::
+   (Num j, LLVM.IsConst j, LLVM.IsInteger j,
+    LLVM.CmpRet j, LLVM.CmpResult j ~ Bool,
+    MultiValue.Additive i, MultiValue.IntegerConstant i) =>
+   LLVM.Value j -> MultiValue.T i -> Iter.T r (MultiValue.T i)
+iteratorStart n start = Iter.take n $ Iter.iterate MultiValue.inc start
+
+
+instance C Word32 where
+   type Index Word32 = Word32
+   intersectCode = MultiValue.min
+   sizeCode (MultiValue.Cons n) = LLVM.ext n
+   size = fromIntegral
+   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
+      Monad.lift2 (,) (LLVM.ext n) (LLVM.ext i)
+   iterator = iteratorPrimitive
+
+instance Sequence Word32 where
+   sequenceShapeFromIndex = return
+
+instance C Word64 where
+   type Index Word64 = Word64
+   intersectCode = MultiValue.min
+   sizeCode (MultiValue.Cons n) = return n
+   size = fromIntegral
+   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) = return (n, i)
+   iterator = iteratorPrimitive
+
+instance Sequence Word64 where
+   sequenceShapeFromIndex = return
+
+
+{- |
+Array dimensions and indexes cannot be negative,
+but computations in indices may temporarily yield negative values
+or we want to add negative values to indices.
+
+Maybe we should better have type Index Word64 = Int64?
+-}
+unsigned8 :: LLVM.Value Int8 -> LLVM.CodeGenFunction r (LLVM.Value Word8)
+unsigned8 = LLVM.bitcast
+
+instance C Int8 where
+   type Index Int8 = Int8
+   intersectCode = MultiValue.min
+   sizeCode (MultiValue.Cons n) = LLVM.ext =<< unsigned8 n
+   size = fromIntegral
+   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
+      Monad.lift2 (,) (LLVM.ext =<< unsigned8 n) (LLVM.ext =<< unsigned8 i)
+   iterator = iteratorPrimitive
+
+instance Sequence Int8 where
+   sequenceShapeFromIndex = return
+
+unsigned16 :: LLVM.Value Int16 -> LLVM.CodeGenFunction r (LLVM.Value Word16)
+unsigned16 = LLVM.bitcast
+
+instance C Int16 where
+   type Index Int16 = Int16
+   intersectCode = MultiValue.min
+   sizeCode (MultiValue.Cons n) = LLVM.ext =<< unsigned16 n
+   size = fromIntegral
+   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
+      Monad.lift2 (,) (LLVM.ext =<< unsigned16 n) (LLVM.ext =<< unsigned16 i)
+   iterator = iteratorPrimitive
+
+instance Sequence Int16 where
+   sequenceShapeFromIndex = return
+
+instance C Int32 where
+   type Index Int32 = Int32
+   intersectCode = MultiValue.min
+   sizeCode (MultiValue.Cons n) = LLVM.zext n
+   size = fromIntegral
+   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
+      Monad.lift2 (,) (LLVM.zext n) (LLVM.zext i)
+   iterator = iteratorPrimitive
+
+instance Sequence Int32 where
+   sequenceShapeFromIndex = return
+
+instance C Int64 where
+   type Index Int64 = Int64
+   intersectCode = MultiValue.min
+   sizeCode (MultiValue.Cons n) = LLVM.bitcast n
+   size = fromIntegral
+   flattenIndexRec (MultiValue.Cons n) (MultiValue.Cons i) =
+      Monad.lift2 (,) (LLVM.bitcast n) (LLVM.bitcast i)
+   iterator = iteratorPrimitive
+
+instance Sequence Int64 where
+   sequenceShapeFromIndex = return
+
+
+{- |
+'Range' denotes an inclusive range like
+those of the Haskell 98 standard @Array@ type from the @array@ package.
+E.g. the shape type @(Range Int32, Range Int64)@
+is equivalent to the ix type @(Int32, Int64)@ for @Array@s.
+-}
+data Range n = Range n n
+
+singletonRange :: n -> Range n
+singletonRange n = Range n n
+
+
+{-# INLINE castToElemPtr #-}
+castToElemPtr :: Ptr (f a) -> Ptr a
+castToElemPtr = castPtr
+
+-- cf. sample-frame:Stereo
+instance Storable n => Storable (Range n) where
+   {-# INLINE sizeOf #-}
+   {-# INLINE alignment #-}
+   {-# INLINE peek #-}
+   {-# INLINE poke #-}
+   sizeOf ~(Range l r) = sizeOf l + mod (- sizeOf l) (alignment r) + sizeOf r
+   alignment ~(Range l _) = alignment l
+   poke p (Range l r) =
+      let q = castToElemPtr p
+      in  poke q l >> pokeElemOff q 1 r
+   peek p =
+      let q = castToElemPtr p
+      in  Monad.lift2 Range (peek q) (peekElemOff q 1)
+
+
+class
+   (MultiValue.Additive n, MultiValue.Real n, MultiValue.IntegerConstant n) =>
+      ToSize n where
+   toSize :: MultiValue.T n -> LLVM.CodeGenFunction r (LLVM.Value Size)
+
+instance ToSize Word32 where toSize (MultiValue.Cons n) = LLVM.ext n
+instance ToSize Word64 where toSize (MultiValue.Cons n) = return n
+instance ToSize Int32 where toSize (MultiValue.Cons n) = LLVM.zext n
+instance ToSize Int64 where toSize (MultiValue.Cons n) = LLVM.bitcast n
+
+rangeSize ::
+   (ToSize n) =>
+   Range (MultiValue.T n) -> LLVM.CodeGenFunction r (LLVM.Value Size)
+rangeSize (Range from to) =
+   toSize =<< MultiValue.inc =<< MultiValue.sub to from
+
+unzipRange :: MultiValue.T (Range n) -> Range (MultiValue.T n)
+unzipRange (MultiValue.Cons (Range from to)) =
+   Range (MultiValue.Cons from) (MultiValue.Cons to)
+
+zipRange :: MultiValue.T n -> MultiValue.T n -> MultiValue.T (Range n)
+zipRange (MultiValue.Cons from) (MultiValue.Cons to) =
+   MultiValue.Cons (Range from to)
+
+instance (MultiValue.C n) => MultiValue.C (Range n) where
+   type Repr f (Range n) = Range (MultiValue.Repr f n)
+   cons (Range from to) = zipRange (MultiValue.cons from) (MultiValue.cons to)
+   undef = MultiValue.compose $ singletonRange MultiValue.undef
+   zero = MultiValue.compose $ singletonRange MultiValue.zero
+   phis bb a =
+      case unzipRange a of
+         Range a0 a1 ->
+            Monad.lift2 zipRange (MultiValue.phis bb a0) (MultiValue.phis bb a1)
+   addPhis bb a b =
+      case (unzipRange a, unzipRange b) of
+         (Range a0 a1, Range b0 b1) ->
+            MultiValue.addPhis bb a0 b0 >>
+            MultiValue.addPhis bb a1 b1
+
+type instance
+   MultiValue.Decomposed f (Range pn) =
+      Range (MultiValue.Decomposed f pn)
+type instance
+   MultiValue.PatternTuple (Range pn) =
+      Range (MultiValue.PatternTuple pn)
+
+instance (MultiValue.Compose n) => MultiValue.Compose (Range n) where
+   type Composed (Range n) = Range (MultiValue.Composed n)
+   compose (Range from to) =
+      zipRange (MultiValue.compose from) (MultiValue.compose to)
+
+instance (MultiValue.Decompose pn) => MultiValue.Decompose (Range pn) where
+   decompose (Range pfrom pto) rng =
+      case unzipRange rng of
+         Range from to ->
+            Range
+               (MultiValue.decompose pfrom from)
+               (MultiValue.decompose pto to)
+
+instance (MultiMem.C n) => MultiMem.C (Range n) where
+   type Struct (Range n) = PairStruct n
+   decompose = fmap (uncurry zipRange) . decomposeGen
+   compose x = case unzipRange x of Range n m -> composeGen n m
+
+instance (Integral n, ToSize n, MultiValue.Comparison n) => C (Range n) where
+   type Index (Range n) = n
+   intersectCode =
+      MultiValue.modifyF2 (singletonRange atom) (singletonRange atom) $
+            \(Range fromN toN) (Range fromM toM) ->
+         Monad.lift2 Range (MultiValue.max fromN fromM) (MultiValue.min toN toM)
+   sizeCode = rangeSize . unzipRange
+   size (Range from to) = fromIntegral $ to-from+1
+   flattenIndexRec rngValue i =
+      case unzipRange rngValue of
+         rng@(Range from _to) ->
+            Monad.lift2 (,) (rangeSize rng) (toSize =<< MultiValue.sub i from)
+   iterator rngValue =
+      case MultiValue.decompose (singletonRange atom) rngValue of
+         Range from to ->
+            IterMV.takeWhile (MultiValue.cmp LLVM.CmpGE to) $
+            Iter.iterate MultiValue.inc from
+
+
+{- |
+'Shifted' denotes a range defined by the start index and the length.
+-}
+data Shifted n = Shifted {shiftedOffset, shiftedSize :: n}
+
+singletonShifted :: n -> Shifted n
+singletonShifted n = Shifted n n
+
+
+-- cf. sample-frame:Stereo
+instance Storable n => Storable (Shifted n) where
+   {-# INLINE sizeOf #-}
+   {-# INLINE alignment #-}
+   {-# INLINE peek #-}
+   {-# INLINE poke #-}
+   sizeOf ~(Shifted l n) = sizeOf l + mod (- sizeOf l) (alignment n) + sizeOf n
+   alignment ~(Shifted l _) = alignment l
+   poke p (Shifted l n) =
+      let q = castToElemPtr p
+      in  poke q l >> pokeElemOff q 1 n
+   peek p =
+      let q = castToElemPtr p
+      in  Monad.lift2 Shifted (peek q) (peekElemOff q 1)
+
+
+unzipShifted :: MultiValue.T (Shifted n) -> Shifted (MultiValue.T n)
+unzipShifted (MultiValue.Cons (Shifted from to)) =
+   Shifted (MultiValue.Cons from) (MultiValue.Cons to)
+
+zipShifted :: MultiValue.T n -> MultiValue.T n -> MultiValue.T (Shifted n)
+zipShifted (MultiValue.Cons from) (MultiValue.Cons to) =
+   MultiValue.Cons (Shifted from to)
+
+instance (MultiValue.C n) => MultiValue.C (Shifted n) where
+   type Repr f (Shifted n) = Shifted (MultiValue.Repr f n)
+   cons (Shifted offset len) =
+      zipShifted (MultiValue.cons offset) (MultiValue.cons len)
+   undef = MultiValue.compose $ singletonShifted MultiValue.undef
+   zero = MultiValue.compose $ singletonShifted MultiValue.zero
+   phis bb a =
+      case unzipShifted a of
+         Shifted a0 a1 ->
+            Monad.lift2 zipShifted
+               (MultiValue.phis bb a0) (MultiValue.phis bb a1)
+   addPhis bb a b =
+      case (unzipShifted a, unzipShifted b) of
+         (Shifted a0 a1, Shifted b0 b1) ->
+            MultiValue.addPhis bb a0 b0 >>
+            MultiValue.addPhis bb a1 b1
+
+type instance
+   MultiValue.Decomposed f (Shifted pn) =
+      Shifted (MultiValue.Decomposed f pn)
+type instance
+   MultiValue.PatternTuple (Shifted pn) =
+      Shifted (MultiValue.PatternTuple pn)
+
+instance (MultiValue.Compose n) => MultiValue.Compose (Shifted n) where
+   type Composed (Shifted n) = Shifted (MultiValue.Composed n)
+   compose (Shifted offset len) =
+      zipShifted (MultiValue.compose offset) (MultiValue.compose len)
+
+instance (MultiValue.Decompose pn) => MultiValue.Decompose (Shifted pn) where
+   decompose (Shifted poffset plen) rng =
+      case unzipShifted rng of
+         Shifted offset len ->
+            Shifted
+               (MultiValue.decompose poffset offset)
+               (MultiValue.decompose plen len)
+
+instance (MultiMem.C n) => MultiMem.C (Shifted n) where
+   type Struct (Shifted n) = PairStruct n
+   decompose = fmap (uncurry zipShifted) . decomposeGen
+   compose x = case unzipShifted x of Shifted n m -> composeGen n m
+
+type PairStruct n = LLVM.Struct (MultiMem.Struct n, (MultiMem.Struct n, ()))
+
+decomposeGen ::
+   (MultiMem.C n) =>
+   LLVM.Value (PairStruct n) ->
+   LLVM.CodeGenFunction r (MultiValue.T n, MultiValue.T n)
+decomposeGen nm =
+   Monad.lift2 (,)
+      (MultiMem.decompose =<< LLVM.extractvalue nm TypeNum.d0)
+      (MultiMem.decompose =<< LLVM.extractvalue nm TypeNum.d1)
+
+composeGen ::
+   (MultiMem.C n) =>
+   MultiValue.T n -> MultiValue.T n ->
+   LLVM.CodeGenFunction r (LLVM.Value (PairStruct n))
+composeGen n m = do
+   sn <- MultiMem.compose n
+   sm <- MultiMem.compose m
+   rn <- LLVM.insertvalue (LLVM.value LLVM.undef) sn TypeNum.d0
+   LLVM.insertvalue rn sm TypeNum.d1
+
+
+instance (Integral n, ToSize n, MultiValue.Comparison n) => C (Shifted n) where
+   type Index (Shifted n) = n
+   intersectCode =
+      MultiValue.modifyF2 (singletonShifted atom) (singletonShifted atom) $
+            \(Shifted offsetN lenN) (Shifted offsetM lenM) -> do
+         offset <- MultiValue.max offsetN offsetM
+         endN <- MultiValue.add offsetN lenN
+         endM <- MultiValue.add offsetM lenM
+         end <- MultiValue.min endN endM
+         Shifted offset <$> MultiValue.sub end offset
+   sizeCode = toSize . shiftedSize . unzipShifted
+   size (Shifted _offset len) = fromIntegral len
+   flattenIndexRec shapeValue i =
+      case unzipShifted shapeValue of
+         Shifted offset len ->
+            Monad.lift2 (,) (toSize len) (toSize =<< MultiValue.sub i offset)
+   iterator rngValue =
+      case MultiValue.decompose (singletonShifted atom) rngValue of
+         Shifted from len ->
+            IterMV.take len $ Iter.iterate MultiValue.inc from
+
+
+
+instance (C n, C m) => C (n,m) where
+   type Index (n,m) = (Index n, Index m)
+   intersectCode a b =
+      case (MultiValue.unzip a, MultiValue.unzip b) of
+         ((an,am), (bn,bm)) ->
+            Monad.lift2 MultiValue.zip
+               (intersectCode an bn)
+               (intersectCode am bm)
+   sizeCode nm =
+      case MultiValue.unzip nm of
+         (n,m) -> Monad.liftJoin2 A.mul (sizeCode n) (sizeCode m)
+   size (n,m) = size n * size m
+   flattenIndexRec nm ij =
+      case (MultiValue.unzip nm, MultiValue.unzip ij) of
+         ((n,m), (i,j)) -> do
+            (ns, il) <- flattenIndexRec n i
+            (ms, jl) <- flattenIndexRec m j
+            Monad.lift2 (,)
+               (A.mul ns ms)
+               (A.add jl =<< A.mul ms il)
+   iterator nm =
+      case MultiValue.unzip nm of
+         (n,m) ->
+            uncurry MultiValue.zip <$>
+            Iter.cartesian (iterator n) (iterator m)
+   loop code nm =
+      case MultiValue.unzip nm of
+         (n,m) -> loop (\i -> loop (\j -> code (MultiValue.zip i j)) m) n
+
+instance (C n, C m, C l) => C (n,m,l) where
+   type Index (n,m,l) = (Index n, Index m, Index l)
+   intersectCode a b =
+      case (MultiValue.unzip3 a, MultiValue.unzip3 b) of
+         ((ai,aj,ak), (bi,bj,bk)) ->
+            Monad.lift3 MultiValue.zip3
+               (intersectCode ai bi)
+               (intersectCode aj bj)
+               (intersectCode ak bk)
+   sizeCode nml =
+      case MultiValue.unzip3 nml of
+         (n,m,l) ->
+            Monad.liftJoin2 A.mul (sizeCode n) $
+            Monad.liftJoin2 A.mul (sizeCode m) (sizeCode l)
+   size (n,m,l) = size n * size m * size l
+   flattenIndexRec nml ijk =
+      case (MultiValue.unzip3 nml, MultiValue.unzip3 ijk) of
+         ((n,m,l), (i,j,k)) -> do
+            (ns, il) <- flattenIndexRec n i
+            (ms, jl) <- flattenIndexRec m j
+            x0 <- A.add jl =<< A.mul ms il
+            (ls, kl) <- flattenIndexRec l k
+            x1 <- A.add kl =<< A.mul ls x0
+            sz <- A.mul ns =<< A.mul ms ls
+            return (sz, x1)
+   iterator nml =
+      case MultiValue.unzip3 nml of
+         (n,m,l) ->
+            fmap (\(a,(b,c)) -> MultiValue.zip3 a b c) $
+            Iter.cartesian (iterator n) $
+            Iter.cartesian (iterator m) (iterator l)
+   loop code nml =
+      case MultiValue.unzip3 nml of
+         (n,m,l) ->
+            loop (\i -> loop (\j -> loop (\k ->
+               code (MultiValue.zip3 i j k))
+            l) m) n
diff --git a/src/Data/Array/Knead/Simple/Fold.hs b/src/Data/Array/Knead/Simple/Fold.hs
--- a/src/Data/Array/Knead/Simple/Fold.hs
+++ b/src/Data/Array/Knead/Simple/Fold.hs
@@ -20,12 +20,12 @@
 import qualified Data.Array.Knead.Simple.Private as Core
 import Data.Array.Knead.Simple.Private (Array(Array), Code, Val, )
 
-import qualified Data.Array.Knead.Index.Linear as Linear
-import qualified Data.Array.Knead.Index.Linear.Int as IndexInt
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Cubic as Linear
+import qualified Data.Array.Knead.Shape.Cubic.Int as IndexInt
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp, unExp, )
-import Data.Array.Knead.Index.Linear ((#:.), (:.)((:.)), )
+import Data.Array.Knead.Shape.Cubic ((#:.), (:.)((:.)), )
 
 import qualified LLVM.Extra.Multi.Value as MultiValue
 import LLVM.Extra.Multi.Value (atom, )
@@ -75,8 +75,7 @@
    (Val (Linear.Index sh) -> Code r a)
 fold1CodeLinear f nc code ix =
    Core.fold1Code f (IndexInt.decons nc)
-      (\jx j -> code (jx #:. IndexInt.cons j))
-      ix
+      (\j -> code (ix #:. IndexInt.cons j))
 
 fold ::
    (MultiValue.C a) =>
diff --git a/src/Data/Array/Knead/Simple/Physical.hs b/src/Data/Array/Knead/Simple/Physical.hs
--- a/src/Data/Array/Knead/Simple/Physical.hs
+++ b/src/Data/Array/Knead/Simple/Physical.hs
@@ -9,7 +9,7 @@
    with,
    render,
    scanl1,
-   mapAccumL,
+   mapAccumLSimple,
    scatter,
    scatterMaybe,
    permute,
@@ -17,7 +17,7 @@
 
 import qualified Data.Array.Knead.Simple.PhysicalPrivate as Priv
 import qualified Data.Array.Knead.Simple.Private as Sym
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import qualified Data.Array.Knead.Code as Code
 import Data.Array.Knead.Expression (Exp, unExp, )
@@ -39,7 +39,6 @@
 
 import Control.Monad.HT (void, )
 import Control.Applicative (liftA2, (<$>), )
-import Data.Word (Word32, )
 
 import Prelude hiding (scanl1, )
 
@@ -101,7 +100,7 @@
 type Importer f = FunPtr f -> f
 
 foreign import ccall safe "dynamic" callShaper ::
-   Importer (Ptr sh -> IO Word32)
+   Importer (Ptr sh -> IO Shape.Size)
 
 foreign import ccall safe "dynamic" callRenderer ::
    Importer (Ptr sh -> Ptr am -> IO ())
@@ -147,18 +146,11 @@
 
 scanl1 ::
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
-    Storable a, MultiValueMemory.C a) =>
-   (Exp a -> Exp a -> Exp a) ->
-   Sym.Array (sh, Word32) a -> IO (Array (sh, Word32) a)
-scanl1 = scanl1Gen
-
-scanl1Gen ::
-   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Shape.C n, Storable n, MultiValueMemory.C n,
     Storable a, MultiValueMemory.C a) =>
    (Exp a -> Exp a -> Exp a) ->
    Sym.Array (sh, n) a -> IO (Array (sh, n) a)
-scanl1Gen f (Sym.Array esh code) =
+scanl1 f (Sym.Array esh code) =
    materialize "scanl1" esh $ \sptr ptr -> do
       (sh, n) <- MultiValue.unzip <$> Shape.load esh sptr
       let step ix ptrStart =
@@ -172,7 +164,7 @@
                 return (ptr1, Maybe.just acc1)
       void $ Shape.loop step sh ptr
 
-mapAccumL ::
+mapAccumLSimple ::
    (Shape.C sh, Storable sh, MultiValueMemory.C sh,
     Shape.C n, Storable n, MultiValueMemory.C n,
     MultiValue.C acc,
@@ -180,9 +172,9 @@
     Storable y, MultiValueMemory.C y) =>
    (Exp acc -> Exp x -> Exp (acc,y)) ->
    Sym.Array sh acc -> Sym.Array (sh, n) x -> IO (Array (sh, n) y)
-mapAccumL f arrInit arrData =
-   materialize "mapAccumL" (Sym.shape arrData) $
-      Priv.mapAccumL f arrInit arrData
+mapAccumLSimple f arrInit arrData =
+   materialize "mapAccumLSimple" (Sym.shape arrData) $
+      Priv.mapAccumLSimple f arrInit arrData
 
 scatterMaybe ::
    (Shape.C sh0, Shape.Index sh0 ~ ix0,
diff --git a/src/Data/Array/Knead/Simple/PhysicalPrivate.hs b/src/Data/Array/Knead/Simple/PhysicalPrivate.hs
--- a/src/Data/Array/Knead/Simple/PhysicalPrivate.hs
+++ b/src/Data/Array/Knead/Simple/PhysicalPrivate.hs
@@ -2,7 +2,7 @@
 module Data.Array.Knead.Simple.PhysicalPrivate where
 
 import qualified Data.Array.Knead.Simple.Private as Sym
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp, unExp)
 import Data.Array.Knead.Code (getElementPtr)
@@ -15,9 +15,9 @@
 
 import qualified LLVM.Core as LLVM
 
-import Foreign.Storable (Storable, )
 import Foreign.Ptr (Ptr, )
 
+import qualified Control.Applicative.HT as App
 import Control.Monad.HT (void, )
 import Control.Applicative ((<$>), )
 
@@ -25,62 +25,177 @@
 
 
 
-mapAccumL ::
-   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
-    Shape.C n, Storable n, MultiValueMemory.C n,
+writeArray ::
+   (Shape.C sh, Memory.C a) =>
+   MultiValue.T sh ->
+   (MultiValue.T (Shape.Index sh) -> LLVM.CodeGenFunction r a) ->
+   LLVM.Value (Ptr (Memory.Struct a)) ->
+   LLVM.CodeGenFunction r (LLVM.Value (Ptr (Memory.Struct a)))
+writeArray sh code ptr = do
+   let clear ix p = do
+         flip Memory.store p =<< code ix
+         A.advanceArrayElementPtr p
+   Shape.loop clear sh ptr
+
+
+mapAccumLLoop ::
+   (MultiValue.C acc,
+    MultiValueMemory.C b, MultiValueMemory.Struct b ~ bm,
+    Shape.C sh, Shape.Index sh ~ ix) =>
+   (MultiValue.T ix -> LLVM.CodeGenFunction r (MultiValue.T a)) ->
+   (Exp acc -> Exp a -> Exp (acc, b)) ->
+   MultiValue.T sh ->
+   LLVM.Value (Ptr bm) -> MultiValue.T acc ->
+   LLVM.CodeGenFunction r (LLVM.Value (Ptr bm), MultiValue.T acc)
+mapAccumLLoop code f n yPtr accInit = do
+   let step k0 (ptr0, acc0) = do
+         x <- code k0
+         (acc1,y) <- MultiValue.unzip <$> Expr.unliftM2 f acc0 x
+         Memory.store y ptr0
+         ptr1 <- A.advanceArrayElementPtr ptr0
+         return (ptr1, acc1)
+   Shape.loop step n (yPtr, accInit)
+
+mapAccumLSimple ::
+   (Shape.C sh, MultiValueMemory.C sh,
+    Shape.C n, MultiValueMemory.C n,
     MultiValue.C acc,
-    Storable x, MultiValueMemory.C x,
-    Storable y, MultiValueMemory.C y) =>
+    MultiValueMemory.C x,
+    MultiValueMemory.C y) =>
    (Exp acc -> Exp x -> Exp (acc,y)) ->
    Sym.Array sh acc -> Sym.Array (sh, n) x ->
    LLVM.Value (Ptr (MultiValueMemory.Struct (sh,n))) ->
    LLVM.Value (Ptr (MultiValueMemory.Struct y)) ->
    LLVM.CodeGenFunction r ()
-mapAccumL f (Sym.Array _ initCode) (Sym.Array esh code) sptr ptr = do
+mapAccumLSimple f (Sym.Array _ initCode) (Sym.Array esh code) sptr ptr = do
    (sh, n) <- MultiValue.unzip <$> Shape.load esh sptr
    let step ix ptrStart = do
          accInit <- initCode ix
-         fmap fst $
-          (\body -> Shape.loop body n (ptrStart, accInit)) $
-                \k0 (ptr0, acc0) -> do
-             x <- code $ MultiValue.zip ix k0
-             (acc1,y) <- MultiValue.unzip <$> Expr.unliftM2 f acc0 x
-             Memory.store y ptr0
-             ptr1 <- A.advanceArrayElementPtr ptr0
-             return (ptr1, acc1)
+         fst <$> mapAccumLLoop (code . MultiValue.zip ix) f n ptrStart accInit
    void $ Shape.loop step sh ptr
 
+mapAccumLSequence ::
+   (Shape.C n, MultiValueMemory.C n,
+    MultiValue.C acc, MultiValueMemory.C final,
+    MultiValueMemory.C x,
+    MultiValueMemory.C y) =>
+   (Exp acc -> Exp x -> Exp (acc,y)) ->
+   (Exp acc -> Exp final) ->
+   Exp acc -> Sym.Array n x ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct final)) ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct n)) ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct y)) ->
+   LLVM.CodeGenFunction r ()
+mapAccumLSequence f final initExp (Sym.Array esh code) accPtr sptr yPtr = do
+   n <- Shape.load esh sptr
+   accInit <- unExp initExp
+   accExit <- snd <$> mapAccumLLoop code f n yPtr accInit
+   flip Memory.store accPtr =<< Expr.unliftM1 final accExit
+
+mapAccumL ::
+   (Shape.C sh, MultiValueMemory.C sh,
+    Shape.C n, MultiValueMemory.C n,
+    MultiValue.C acc, MultiValueMemory.C final,
+    MultiValueMemory.C x,
+    MultiValueMemory.C y) =>
+   (Exp acc -> Exp x -> Exp (acc,y)) ->
+   (Exp acc -> Exp final) ->
+   Sym.Array sh acc -> Sym.Array (sh, n) x ->
+   (LLVM.Value (Ptr (MultiValueMemory.Struct sh)),
+    LLVM.Value (Ptr (MultiValueMemory.Struct final))) ->
+   (LLVM.Value (Ptr (MultiValueMemory.Struct (sh,n))),
+    LLVM.Value (Ptr (MultiValueMemory.Struct y))) ->
+   LLVM.CodeGenFunction r ()
+mapAccumL f final (Sym.Array _ initCode) (Sym.Array esh code)
+      (_, accPtr) (sptr, yPtr) = do
+   (sh, n) <- MultiValue.unzip <$> Shape.load esh sptr
+   let step ix (accPtr0, yPtrStart) = do
+         accInit <- initCode ix
+         (ptrStop, accExit) <-
+            mapAccumLLoop (code . MultiValue.zip ix) f n yPtrStart accInit
+         flip Memory.store accPtr0 =<< Expr.unliftM1 final accExit
+         accPtr1 <- A.advanceArrayElementPtr accPtr0
+         return (accPtr1, ptrStop)
+   void $ Shape.loop step sh (accPtr,yPtr)
+
 foldOuterL ::
-   (Shape.C sh, Storable sh, MultiValueMemory.C sh,
-    Shape.C n, Storable n, MultiValueMemory.C n,
+   (Shape.C sh, MultiValueMemory.C sh,
+    Shape.C n, MultiValueMemory.C n,
     MultiValueMemory.C a) =>
    (Exp a -> Exp b -> Exp a) ->
    Sym.Array sh a -> Sym.Array (n,sh) b ->
    LLVM.Value (Ptr (MultiValueMemory.Struct sh)) ->
    LLVM.Value (Ptr (MultiValueMemory.Struct a)) ->
    LLVM.CodeGenFunction r ()
-foldOuterL f (Sym.Array _ initCode) (Sym.Array esh code) _sptr ptr = do
-   -- (n,sh) <- MultiValue.unzip <$> Shape.load esh sptr
-   (n,sh) <- MultiValue.unzip <$> unExp esh
-   let fillInit ix ptr0 = do
-         a <- initCode ix
-         Memory.store a ptr0
-         A.advanceArrayElementPtr ptr0
-   void $ Shape.loop fillInit sh ptr
+foldOuterL f (Sym.Array _ initCode) (Sym.Array esh code) sptr ptr = do
+   sh <- Shape.load (Expr.snd esh) sptr
+   n <- MultiValue.fst <$> unExp esh
+   void $ writeArray sh initCode ptr
 
    let step k ix ptr0 = do
-       b <- code $ MultiValue.zip k ix
-       a0 <- Memory.load ptr0
-       a1 <- Expr.unliftM2 f a0 b
-       Memory.store a1 ptr0
-       A.advanceArrayElementPtr ptr0
+         b <- code $ MultiValue.zip k ix
+         a0 <- Memory.load ptr0
+         a1 <- Expr.unliftM2 f a0 b
+         Memory.store a1 ptr0
+         A.advanceArrayElementPtr ptr0
    void $ Shape.loop (\k () -> void $ Shape.loop (step k) sh ptr) n ()
 
+{- |
+We need a scalar Shape type @n@.
+Scalar Shape types could be distinguished from other Shape types
+by the fact that you can convert any Index into a Shape.
+-}
+mapFilter ::
+   (Shape.Sequence n, MultiValueMemory.C n,
+    MultiValueMemory.C b) =>
+   (Exp a -> Exp b) ->
+   (Exp a -> Exp Bool) ->
+   Sym.Array n a ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct n)) ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct b)) ->
+   LLVM.CodeGenFunction r (MultiValue.T n)
+mapFilter f p (Sym.Array esh code) sptr ptr = do
+   n <- Shape.load esh sptr
+   let step ix (dstPtr,dstIx) = do
+         a <- code ix
+         MultiValue.Cons c <- Expr.unliftM1 p a
+         C.ifThen c (dstPtr,dstIx)
+            (do
+               flip Memory.store dstPtr =<< Expr.unliftM1 f a
+               App.lift2 (,)
+                  (A.advanceArrayElementPtr dstPtr)
+                  (MultiValue.inc dstIx))
+   Shape.sequenceShapeFromIndex . snd
+      =<< Shape.loop step n (ptr, MultiValue.zero)
+
+filterOuter ::
+   (Shape.Sequence n, MultiValueMemory.C n,
+    Shape.C sh, MultiValueMemory.C sh,
+    MultiValueMemory.C a) =>
+   Sym.Array n Bool ->
+   Sym.Array (n,sh) a ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct (n,sh))) ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct a)) ->
+   LLVM.CodeGenFunction r (MultiValue.T (n,sh))
+filterOuter (Sym.Array _eish selectCode) (Sym.Array esh code) sptr ptr = do
+   (n,sh) <- MultiValue.unzip <$> Shape.load esh sptr
+   let step k (dstPtr0,dstK) = do
+         MultiValue.Cons c <- selectCode k
+         C.ifThen c (dstPtr0,dstK)
+            (do
+               dstPtr1 <- writeArray sh (code . MultiValue.zip k) dstPtr0
+               (,) dstPtr1 <$> MultiValue.inc dstK)
+   finalN <-
+      Shape.sequenceShapeFromIndex . snd
+         =<< Shape.loop step n (ptr, MultiValue.zero)
+   return $ MultiValue.zip finalN sh
+
+
 scatterMaybe ::
    (Shape.C sh0, Shape.Index sh0 ~ ix0,
     Shape.C sh1, Shape.Index sh1 ~ ix1,
-    Storable sh1, MultiValueMemory.C sh1,
-    Storable a, MultiValueMemory.C a) =>
+    MultiValueMemory.C sh1,
+    MultiValueMemory.C a) =>
    (Exp a -> Exp a -> Exp a) ->
    Sym.Array sh1 a -> Sym.Array sh0 (Maybe (ix1, a)) ->
    LLVM.Value (Ptr (MultiValueMemory.Struct sh1)) ->
@@ -89,11 +204,8 @@
 scatterMaybe accum (Sym.Array esh codeInit) (Sym.Array eish codeMap)
       sptr ptr = do
 
-   let clear ix p = do
-         flip Memory.store p =<< codeInit ix
-         A.advanceArrayElementPtr p
    sh <- Shape.load esh sptr
-   void $ Shape.loop clear sh ptr
+   void $ writeArray sh codeInit ptr
 
    ish <- unExp eish
    let fill ix () = do
@@ -109,8 +221,8 @@
 scatter ::
    (Shape.C sh0, Shape.Index sh0 ~ ix0,
     Shape.C sh1, Shape.Index sh1 ~ ix1,
-    Storable sh1, MultiValueMemory.C sh1,
-    Storable a, MultiValueMemory.C a) =>
+    MultiValueMemory.C sh1,
+    MultiValueMemory.C a) =>
    (Exp a -> Exp a -> Exp a) ->
    Sym.Array sh1 a ->
    Sym.Array sh0 (Shape.Index sh1, a) ->
@@ -118,17 +230,32 @@
    LLVM.Value (Ptr (MultiValueMemory.Struct a)) ->
    LLVM.CodeGenFunction r ()
 scatter accum (Sym.Array esh codeInit) (Sym.Array eish codeMap) sptr ptr = do
-   let clear ix p = do
-         flip Memory.store p =<< codeInit ix
-         A.advanceArrayElementPtr p
    sh <- Shape.load esh sptr
-   void $ Shape.loop clear sh ptr
+   void $ writeArray sh codeInit ptr
 
    ish <- unExp eish
    let fill ix () = do
-         (jx, a) <- fmap MultiValue.unzip $ codeMap ix
+         (jx, a) <- MultiValue.unzip <$> codeMap ix
          p <- getElementPtr sh ptr jx
          flip Memory.store p
             =<< Expr.unliftM2 (flip accum) a
             =<< Memory.load p
    Shape.loop fill ish ()
+
+addDimension ::
+   (Shape.C n, MultiValueMemory.C n, Shape.Index n ~ k,
+    Shape.C sh, MultiValueMemory.C sh,
+    MultiValueMemory.C b) =>
+   Exp n ->
+   (Exp k -> Exp a -> Exp b) ->
+   Sym.Array sh a ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct (sh,n))) ->
+   LLVM.Value (Ptr (MultiValueMemory.Struct b)) ->
+   LLVM.CodeGenFunction r ()
+addDimension en select (Sym.Array esh code) sptr ptr = do
+   (sh,n) <- MultiValue.unzip <$> Shape.load (Expr.zip esh en) sptr
+
+   let fill ix ptr0 = do
+         a <- code ix
+         writeArray n (\k -> Expr.unliftM2 select k a) ptr0
+   void $ Shape.loop fill sh ptr
diff --git a/src/Data/Array/Knead/Simple/Private.hs b/src/Data/Array/Knead/Simple/Private.hs
--- a/src/Data/Array/Knead/Simple/Private.hs
+++ b/src/Data/Array/Knead/Simple/Private.hs
@@ -2,17 +2,17 @@
 {-# LANGUAGE TypeFamilies #-}
 module Data.Array.Knead.Simple.Private where
 
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp(Exp), )
 
 import qualified LLVM.Extra.Multi.Value as MultiValue
-import qualified LLVM.Extra.Monad as Monad
+import qualified LLVM.Extra.Iterator as Iter
 import qualified LLVM.Extra.Maybe as Maybe
-import qualified LLVM.Extra.Control as C
 import qualified LLVM.Core as LLVM
 
 import qualified Control.Category as Cat
+import qualified Control.Monad.HT as Monad
 import Control.Monad ((<=<), )
 
 import Prelude hiding (id, map, zipWith, replicate, )
@@ -86,15 +86,15 @@
    lift2 $ \(Array _sha codeA) (Array _shb codeB) ->
       Array sh
          (\ix ->
-            Monad.liftR2 (Expr.unliftM2 f)
+            Monad.liftJoin2 (Expr.unliftM2 f)
                (codeA =<< Expr.unliftM1 projectIndex0 ix)
                (codeB =<< Expr.unliftM1 projectIndex1 ix))
 
 
 id ::
-   (Shape.C sh, Shape.Index sh ~ ix) =>
-   Exp sh -> Array sh ix
-id sh = Array sh return
+   (C array, Shape.C sh, Shape.Index sh ~ ix) =>
+   Exp sh -> array sh ix
+id sh = lift0 $ Array sh return
 
 map ::
    (C array, Shape.C sh) =>
@@ -114,17 +114,17 @@
 
 
 fold1Code ::
-   (Shape.C sh1, Shape.Index sh1 ~ ix1, MultiValue.C a) =>
+   (Shape.C sh, Shape.Index sh ~ ix, MultiValue.C a) =>
    (Exp a -> Exp a -> Exp a) ->
-   Exp sh1 ->
-   (Val ix0 -> Val ix1 -> Code r a) ->
-   (Val ix0 -> Code r a)
-fold1Code f (Exp nc) code ix = do
+   Exp sh ->
+   (Val ix -> Code r a) ->
+   Code r a
+fold1Code f (Exp nc) code = do
    n <- nc
    fmap Maybe.fromJust $
       Shape.loop
          (\i0 macc0 -> do
-            a <- code ix i0
+            a <- code i0
             acc1 <- Maybe.run macc0 (return a) (flip (Expr.unliftM2 f) a)
             return $ Maybe.just acc1)
          n Maybe.nothing
@@ -136,25 +136,16 @@
 fold1 f =
    lift1 $ \(Array shs code) ->
       case Expr.unzip shs of
-         (sh, s) -> Array sh $ fold1Code f s $ MultiValue.curry code
+         (sh, s) -> Array sh $ fold1Code f s . MultiValue.curry code
 
 
 fold1All ::
    (Shape.C sh, MultiValue.C a) =>
    (Exp a -> Exp a -> Exp a) ->
-   Array sh a -> Array () a
-fold1All f (Array esh code) =
-   fold1 f $
-   Array
-      (Expr.lift1 (MultiValue.zip (MultiValue.Cons ())) esh)
-      (code . MultiValue.snd)
+   Array sh a -> Exp a
+fold1All f (Array sh code) = Exp (fold1Code f sh code)
 
 
-{-
-'Shape.loop' does not support an early exit.
-I wished the LLVM optimizer would rewrite the loop accordingly.
-Unfortunately, it does not.
--}
 findAllCode ::
    (Shape.C sh, Shape.Index sh ~ ix, MultiValue.C a) =>
    (Exp a -> Exp Bool) ->
@@ -164,23 +155,20 @@
 findAllCode p (Exp sh) code = do
    n <- sh
    finalFound <-
-      Shape.loop
-         (\i found ->
-            C.ifThenElse (Maybe.isJust found)
-               (return found)
-               (do
-                  a <- code i
-                  MultiValue.Cons b <- Expr.unliftM1 p a
-                  return $ Maybe.fromBool b a))
-         n Maybe.nothing
+      Iter.mapWhileState_
+         (\a _found -> do
+            MultiValue.Cons b <- Expr.unliftM1 p a
+            notb <- LLVM.inv b
+            return (notb, Maybe.fromBool b a))
+         (Iter.mapM code $ Shape.iterator n)
+         Maybe.nothing
    Maybe.run finalFound
       (return MultiValue.nothing)
       (return . MultiValue.just)
 
 {- |
 In principle this can be implemented using fold1All
-but this one should have short-cut semantics.
-Currently it has not! :-(
+but it has a short-cut semantics.
 @All@ means that it scans all dimensions
 but it does not mean that it finds all occurrences.
 If you want to get the index of the found element,
diff --git a/src/Data/Array/Knead/Simple/ShapeDependent.hs b/src/Data/Array/Knead/Simple/ShapeDependent.hs
--- a/src/Data/Array/Knead/Simple/ShapeDependent.hs
+++ b/src/Data/Array/Knead/Simple/ShapeDependent.hs
@@ -4,12 +4,11 @@
 import qualified Data.Array.Knead.Simple.Private as Core
 import Data.Array.Knead.Simple.Private (Array(Array), )
 
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp, )
 
-import qualified LLVM.Extra.Monad as Monad
-
+import qualified Control.Monad.HT as Monad
 import Control.Monad ((<=<), )
 
 
@@ -62,7 +61,7 @@
    Core.lift2 $ \(Array sha codeA) (Array shb codeB) ->
       Array (combineShape sha shb)
          (\ix ->
-            Monad.liftR2 (Expr.unliftM2 f)
+            Monad.liftJoin2 (Expr.unliftM2 f)
                (codeA =<< Expr.unliftM1 projectIndex0 ix)
                (codeB =<< Expr.unliftM1 projectIndex1 ix))
 
diff --git a/src/Data/Array/Knead/Simple/Slice.hs b/src/Data/Array/Knead/Simple/Slice.hs
--- a/src/Data/Array/Knead/Simple/Slice.hs
+++ b/src/Data/Array/Knead/Simple/Slice.hs
@@ -69,10 +69,10 @@
 import qualified Data.Array.Knead.Simple.ShapeDependent as ShapeDep
 import qualified Data.Array.Knead.Simple.Private as Core
 
-import qualified Data.Array.Knead.Index.Linear as Linear
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Cubic as Linear
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
-import Data.Array.Knead.Index.Linear ((#:.), (:.)((:.)), )
+import Data.Array.Knead.Shape.Cubic ((#:.), (:.)((:.)), )
 import Data.Array.Knead.Expression (Exp, )
 
 import qualified LLVM.Extra.Multi.Value as MultiValue
diff --git a/src/Data/Array/Knead/Simple/Symbolic.hs b/src/Data/Array/Knead/Simple/Symbolic.hs
--- a/src/Data/Array/Knead/Simple/Symbolic.hs
+++ b/src/Data/Array/Knead/Simple/Symbolic.hs
@@ -4,6 +4,7 @@
    Core.Array,
    Core.C(..),
    Exp,
+   fix,
    shape,
    (Core.!),
    Core.the,
@@ -21,24 +22,29 @@
    zip,
    zip3,
    zip4,
-   fold1,
-   fold1All,
+   Core.fold1,
+   Core.fold1All,
    Core.findAll,
    ) where
 
 import qualified Data.Array.Knead.Simple.ShapeDependent as ShapeDep
 import qualified Data.Array.Knead.Simple.Private as Core
-import Data.Array.Knead.Simple.Private (Array, shape, fold1, gather, )
+import Data.Array.Knead.Simple.Private (Array, shape, gather, )
 
-import qualified Data.Array.Knead.Index.Nested.Shape as Shape
+import qualified Data.Array.Knead.Shape.Nested as Shape
 import qualified Data.Array.Knead.Expression as Expr
 import Data.Array.Knead.Expression (Exp, )
 
 import qualified LLVM.Extra.Multi.Value as MultiValue
 
+import Data.Function.HT (Id)
+
 import Prelude hiding (zipWith, zipWith3, zip, zip3, replicate, )
 
 
+fix :: Id (Array sh a)
+fix = id
+
 backpermute ::
    (Shape.C sh0, Shape.Index sh0 ~ ix0,
     Shape.C sh1, Shape.Index sh1 ~ ix1,
@@ -85,16 +91,3 @@
    array sh a -> array sh b -> array sh c -> array sh d ->
    array sh (a,b,c,d)
 zip4 = zipWith4 (Expr.lift4 MultiValue.zip4)
-
-
-fold1All ::
-   (Core.C array, Shape.C sh, Shape.Scalar z, MultiValue.C a) =>
-   (Exp a -> Exp a -> Exp a) ->
-   array sh a -> array z a
-fold1All f =
-   Core.lift1 $ \arr ->
-      fold1 f $
-      backpermute
-         (Expr.lift2 MultiValue.zip Shape.scalar (shape arr))
-         (Expr.lift1 MultiValue.snd)
-         arr