diff --git a/egison.cabal b/egison.cabal
--- a/egison.cabal
+++ b/egison.cabal
@@ -1,5 +1,5 @@
 Name:                egison
-Version:             3.7.10
+Version:             3.7.11
 Synopsis:            Programming language with non-linear pattern-matching against non-free data
 Description:
   An interpreter for Egison, a **pattern-matching-oriented**, purely functional programming language.
@@ -68,7 +68,7 @@
   location: https://github.com/egison/egison.git
   
 Library
-  Build-Depends:   base >= 4.0 && < 5, array, random, containers, unordered-containers, haskeline, transformers, mtl, parsec >= 3.0, directory, ghc, ghc-paths, text, regex-tdfa, process, vector, parallel
+  Build-Depends:   base >= 4.0 && < 5, array, random, containers, unordered-containers, haskeline, transformers, mtl, parsec >= 3.0, directory, ghc, ghc-paths, text, regex-tdfa, process, vector, parallel, split, hashable
   Hs-Source-Dirs:  hs-src
   Exposed-Modules:
                    Language.Egison
diff --git a/hs-src/Interpreter/egison.hs b/hs-src/Interpreter/egison.hs
--- a/hs-src/Interpreter/egison.hs
+++ b/hs-src/Interpreter/egison.hs
@@ -71,7 +71,7 @@
                                                   else evalEgisonTopExprsTestOnly env [LoadFile file]
                                       either print (const $ return ()) result
                                     Options {optTestOnly = False} -> do
-                                      result <- evalEgisonTopExprs env [LoadFile file, Execute (ApplyExpr (VarExpr "main") (CollectionExpr (map (ElementExpr . StringExpr) (map T.pack args))))]
+                                      result <- evalEgisonTopExprs env [LoadFile file, Execute (ApplyExpr (VarExpr $ stringToVar "main") (CollectionExpr (map (ElementExpr . StringExpr) (map T.pack args))))]
                                       either print (const $ return ()) result
 
 data Options = Options {
diff --git a/hs-src/Language/Egison/Core.hs b/hs-src/Language/Egison/Core.hs
--- a/hs-src/Language/Egison/Core.hs
+++ b/hs-src/Language/Egison/Core.hs
@@ -49,6 +49,7 @@
 import Control.Monad.Trans.Maybe
 
 import Data.List (partition)
+import Data.List.Split (split, oneOf)
 import Data.Sequence (Seq, ViewL(..), ViewR(..), (><))
 import qualified Data.Sequence as Sq
 import Data.Ratio
@@ -81,10 +82,10 @@
   forM_ rest $ evalTopExpr env
   return env
  where
-  collectDefs :: [EgisonTopExpr] -> [(String, EgisonExpr)] -> [EgisonTopExpr] -> EgisonM ([(String, EgisonExpr)], [EgisonTopExpr])
+  collectDefs :: [EgisonTopExpr] -> [(Var, EgisonExpr)] -> [EgisonTopExpr] -> EgisonM ([(Var, EgisonExpr)], [EgisonTopExpr])
   collectDefs (expr:exprs) bindings rest =
     case expr of
-      Define name expr -> collectDefs exprs ((show name, expr) : bindings) rest
+      Define name expr -> collectDefs exprs ((name, expr) : bindings) rest
       Load file -> do
         exprs' <- loadLibraryFile file
         collectDefs (exprs' ++ exprs) bindings rest
@@ -102,10 +103,10 @@
   forM_ rest $ evalTopExpr env
   return env
  where
-  collectDefs :: [EgisonTopExpr] -> [(String, EgisonExpr)] -> [EgisonTopExpr] -> EgisonM ([(String, EgisonExpr)], [EgisonTopExpr])
+  collectDefs :: [EgisonTopExpr] -> [(Var, EgisonExpr)] -> [EgisonTopExpr] -> EgisonM ([(Var, EgisonExpr)], [EgisonTopExpr])
   collectDefs (expr:exprs) bindings rest =
     case expr of
-      Define name expr -> collectDefs exprs ((show name, expr) : bindings) rest
+      Define name expr -> collectDefs exprs ((name, expr) : bindings) rest
       Load file -> do
         exprs' <- loadLibraryFile file
         collectDefs (exprs' ++ exprs) bindings rest
@@ -124,10 +125,10 @@
   forM_ rest $ evalTopExpr env
   return env
  where
-  collectDefs :: [EgisonTopExpr] -> [(String, EgisonExpr)] -> [EgisonTopExpr] -> EgisonM ([(String, EgisonExpr)], [EgisonTopExpr])
+  collectDefs :: [EgisonTopExpr] -> [(Var, EgisonExpr)] -> [EgisonTopExpr] -> EgisonM ([(Var, EgisonExpr)], [EgisonTopExpr])
   collectDefs (expr:exprs) bindings rest =
     case expr of
-      Define name expr -> collectDefs exprs ((show name, expr) : bindings) rest
+      Define name expr -> collectDefs exprs ((name, expr) : bindings) rest
       Load _ -> throwError $ Default "No IO support"
       LoadFile _ -> throwError $ Default "No IO support"
       _ -> collectDefs exprs bindings (expr : rest)
@@ -142,8 +143,8 @@
   return $ snd ret
 
 evalTopExpr' :: Env -> EgisonTopExpr -> EgisonM (Maybe String, Env)
-evalTopExpr' env (Define name expr) = recursiveBind env [(show name, expr)] >>= return . ((,) Nothing)
-evalTopExpr' env (Redefine name expr) = recursiveRebind env (show name, expr) >>= return . ((,) Nothing)
+evalTopExpr' env (Define name expr) = recursiveBind env [(name, expr)] >>= return . ((,) Nothing)
+evalTopExpr' env (Redefine name expr) = recursiveRebind env (name, expr) >>= return . ((,) Nothing)
 evalTopExpr' env (Test expr) = do
   val <- evalExprDeep env expr
   return (Just (show val), env)
@@ -168,19 +169,26 @@
     Value (ScalarData s) -> return . Value $ ScalarData $ Div (Plus [Term 1 [(Quote s, 1)]]) (Plus [Term 1 []])
     _ -> throwError $ TypeMismatch "scalar in quote" $ whnf
 
-evalExpr env (QuoteFunctionExpr expr) = do
+evalExpr env (QuoteSymbolExpr expr) = do
   whnf <- evalExpr env expr
   case whnf of
     Value val -> return . Value $ QuotedFunc val
     _ -> throwError $ TypeMismatch "value in quote-function" $ whnf
 
-evalExpr env (VarExpr name) = refVar' env name >>= evalRef
+evalExpr env (VarExpr name) = do
+  x <- refVar' env name >>= evalRef
+  return (case x of
+            Value (ScalarData (Div (Plus [Term 1 [(FunctionData fn argnames args js, 1)]]) p)) -> 
+              case fn of
+                Nothing -> Value $ ScalarData (Div (Plus [Term 1 [(FunctionData (Just $ symbolScalarData "" $ show name) argnames args js, 1)]]) p)
+                Just s -> Value $ ScalarData (Div (Plus [Term 1 [(FunctionData fn argnames args js, 1)]]) p)
+            _ -> x)
  where
-  refVar' :: Env -> String -> EgisonM ObjectRef
-  refVar' env var = maybe (newEvaluatedObjectRef (Value (symbolScalarData "" var))) return
+  refVar' :: Env -> Var -> EgisonM ObjectRef
+  refVar' env var = maybe (newEvaluatedObjectRef (Value (symbolScalarData "" $ show var))) return
                           (refVar env var)
 
-evalExpr env (PartialVarExpr n) = evalExpr env (VarExpr ("::" ++ show n))
+evalExpr env (PartialVarExpr n) = evalExpr env (VarExpr $ stringToVar ("::" ++ show n))
 
 evalExpr _ (InductiveDataExpr name []) = return . Value $ InductiveData name []
 evalExpr env (InductiveDataExpr name exprs) =
@@ -200,7 +208,7 @@
   fromInnerExpr :: InnerExpr -> EgisonM Inner
   fromInnerExpr (ElementExpr expr) = IElement <$> newObjectRef env expr
   fromInnerExpr (SubCollectionExpr expr) = ISubCollection <$> newObjectRef env expr
-
+ 
 evalExpr env (ArrayExpr exprs) = do
   refs' <- mapM (newObjectRef env) exprs
   return . Intermediate . IArray $ Array.listArray (1, toInteger (length exprs)) refs'
@@ -257,18 +265,10 @@
       _ -> throwError $ TypeMismatch "integer or string" $ Value val
   makeHashKey whnf = throwError $ TypeMismatch "integer or string" $ whnf
 
-evalExpr env (UserIndexedExpr expr indices) = do
-  val <- evalExprDeep env expr
-  js <- mapM (\i -> case i of
-                      Userscript n -> evalExprDeep env n >>= return . Userscript) indices
-  case val of
-    (UserIndexedData val' is') -> return $ Value $ UserIndexedData val' (is' ++ js)
-    _ -> return $ Value $ UserIndexedData val js
-
-evalExpr env (IndexedExpr False expr indices) = do
+evalExpr env (IndexedExpr bool expr indices) = do
   tensor <- case expr of
-              (VarExpr var) -> do
-                let mObjRef = refVar env (show (Var var (map f indices)))
+              VarExpr (Var xs is) -> do
+                let mObjRef = refVar env (Var xs $ is ++ (map f indices))
                 case mObjRef of
                   (Just objRef) -> evalRef objRef
                   Nothing -> evalExpr env expr
@@ -277,6 +277,7 @@
                       Superscript n -> evalExprDeep env n >>= return . Superscript
                       Subscript n -> evalExprDeep env n >>= return . Subscript
                       SupSubscript n -> evalExprDeep env n >>= return . SupSubscript
+                      Userscript n -> evalExprDeep env n >>= return . Userscript
               ) indices
   
   ret <- case tensor of
@@ -285,37 +286,58 @@
                              Superscript n -> evalExprDeep env n >>= extractScalar >>= return . Superscript
                              Subscript n -> evalExprDeep env n >>= extractScalar >>= return . Subscript
                              SupSubscript n -> evalExprDeep env n >>= extractScalar >>= return . SupSubscript
+                             Userscript n -> evalExprDeep env n >>= extractScalar >>= return . Userscript
                     ) indices
         return $ Value (ScalarData (Div (Plus [(Term 1 [(Symbol id name js2, 1)])]) (Plus [(Term 1 [])])))
-      (Value (ScalarData _)) -> do
-        return $ tensor
-      (Value (TensorData (Tensor ns xs _))) -> do
-        tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
-      (Intermediate (ITensor (Tensor ns xs _))) -> do
-        tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor
+      (Value (ScalarData (Div (Plus [(Term 1 [(Symbol id name js', 1)])]) (Plus [(Term 1 [])])))) -> do
+        js2 <- mapM (\i -> case i of
+                             Superscript n -> evalExprDeep env n >>= extractScalar >>= return . Superscript
+                             Subscript n -> evalExprDeep env n >>= extractScalar >>= return . Subscript
+                             SupSubscript n -> evalExprDeep env n >>= extractScalar >>= return . SupSubscript
+                             Userscript n -> evalExprDeep env n >>= extractScalar >>= return . Userscript
+                    ) indices
+        return $ Value (ScalarData (Div (Plus [(Term 1 [(Symbol id name (js' ++ js2), 1)])]) (Plus [(Term 1 [])])))
+      (Value (TensorData (Tensor ns xs is))) -> do
+        if bool then tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
+                else tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
+      (Intermediate (ITensor (Tensor ns xs is))) -> do
+        if bool then tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor
+                else tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor
       _ -> do
         js2 <- mapM (\i -> case i of
                              Superscript n -> evalExprDeep env n >>= extractScalar >>= return . Superscript
                              Subscript n -> evalExprDeep env n >>= extractScalar >>= return . Subscript
                              SupSubscript n -> evalExprDeep env n >>= extractScalar >>= return . SupSubscript
+                             Userscript n -> evalExprDeep env n >>= extractScalar >>= return . Userscript
                     ) indices
         refArray tensor (map (\j -> case j of
                                       Superscript k -> ScalarData k
                                       Subscript k -> ScalarData k
                                       SupSubscript k -> ScalarData k
+                                      Userscript k -> ScalarData k
                               ) js2)
-  return ret
+  let ret2 = case expr of
+               (VarExpr var) -> do
+                 case ret of
+                   Value (ScalarData (Div (Plus [Term 1 [(FunctionData fn argnames args js, 1)]]) p)) ->
+                     case fn of
+                       Nothing -> Value $ ScalarData (Div (Plus [Term 1 [(FunctionData (Just $ symbolScalarData "" $ show var ++ concat (map show indices)) argnames args js, 1)]]) p)
+                       Just s -> Value $ ScalarData (Div (Plus [Term 1 [(FunctionData fn argnames args js, 1)]]) p)
+                   _ -> ret
+               _ -> ret
+  return ret2
  where
   f :: Index a -> Index ()
   f (Superscript _) = Superscript ()
   f (Subscript _) = Subscript ()
   f (SupSubscript _) = SupSubscript ()
+  f (Userscript _) = Userscript ()
 
-evalExpr env (SubrefsExpr expr jsExpr) = do
+evalExpr env (SubrefsExpr bool expr jsExpr) = do
   js <- evalExpr env jsExpr >>= collectionToList >>= return . (map Subscript)
   tensor <- case expr of
-              (VarExpr var) -> do
-                let mObjRef = refVar env (show (Var var (take (length js) (repeat (Subscript ())))))
+              VarExpr (Var xs is) -> do
+                let mObjRef = refVar env (Var xs $ is ++ (take (length js) (repeat (Subscript ()))))
                 case mObjRef of
                   (Just objRef) -> evalRef objRef
                   Nothing -> evalExpr env expr
@@ -324,9 +346,11 @@
       (Value (ScalarData _)) -> do
         return $ tensor
       (Value (TensorData (Tensor ns xs is))) -> do
-        tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
+        if bool then tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
+                else tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
       (Intermediate (ITensor (Tensor ns xs is))) -> do
-        tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor
+        if bool then tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor
+                else tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor
       _ -> throwError $ NotImplemented "subrefs"
   return ret
  where
@@ -334,12 +358,13 @@
   f (Superscript _) = Superscript ()
   f (Subscript _) = Subscript ()
   f (SupSubscript _) = SupSubscript ()
+  f (Userscript _) = Userscript ()
 
-evalExpr env (SuprefsExpr expr jsExpr) = do
+evalExpr env (SuprefsExpr bool expr jsExpr) = do
   js <- evalExpr env jsExpr >>= collectionToList >>= return . (map Superscript)
   tensor <- case expr of
-              (VarExpr var) -> do
-                let mObjRef = refVar env (show (Var var (take (length js) (repeat (Superscript ())))))
+              VarExpr (Var xs is) -> do
+                let mObjRef = refVar env (Var xs $ is ++ (take (length js) (repeat (Superscript ()))))
                 case mObjRef of
                   (Just objRef) -> evalRef objRef
                   Nothing -> evalExpr env expr
@@ -348,9 +373,11 @@
       (Value (ScalarData _)) -> do
         return $ tensor
       (Value (TensorData (Tensor ns xs is))) -> do
-        tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
+        if bool then tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
+                else tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor >>= return . Value
       (Intermediate (ITensor (Tensor ns xs is))) -> do
-        tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor
+        if bool then tref js (Tensor ns xs js) >>= toTensor >>= tContract' >>= fromTensor
+                else tref (is ++ js) (Tensor ns xs (is ++ js)) >>= toTensor >>= tContract' >>= fromTensor
       _ -> throwError $ NotImplemented "suprefs"
   return ret
  where
@@ -358,7 +385,17 @@
   f (Superscript _) = Superscript ()
   f (Subscript _) = Subscript ()
   f (SupSubscript _) = SupSubscript ()
+  f (Userscript _) = Userscript ()
 
+evalExpr env (UserrefsExpr bool expr jsExpr) = do
+  val <- evalExprDeep env expr
+  js <- evalExpr env jsExpr >>= collectionToList >>= mapM extractScalar >>= return . (map Userscript)
+  ret <- case val of
+      (ScalarData (Div (Plus [Term 1 [(Symbol id name is, 1)]]) (Plus [Term 1 []]))) -> return $ Value (ScalarData (Div (Plus [Term 1 [(Symbol id name (is ++ js), 1)]]) (Plus [Term 1 []])))
+      (ScalarData (Div (Plus [Term 1 [(FunctionData (Just name) argnames args is, 1)]]) (Plus [Term 1 []]))) -> return $ Value (ScalarData (Div (Plus [Term 1 [(FunctionData (Just name) argnames args (is ++ js), 1)]]) (Plus [Term 1 []])))
+      _ -> throwError $ NotImplemented "user-refs"
+  return ret
+
 evalExpr env (LambdaExpr names expr) = do
   names' <- mapM (\name -> case name of
                              (TensorArg name') -> return name'
@@ -375,6 +412,21 @@
 
 evalExpr env (PatternFunctionExpr names pattern) = return . Value $ PatternFunc env names pattern
 
+evalExpr (Env frame Nothing) (FunctionExpr args) = throwError $ Default "function symbol is not bound to a variable" 
+
+evalExpr env@(Env frame (Just name)) (FunctionExpr args) = do
+  args' <- mapM (\arg -> evalExprDeep env arg) args
+  return . Value $ ScalarData (Div (Plus [Term 1 [(FunctionData (Just $ symbolScalarData "" $ show name) (map (\x -> symbolScalarData "" $ show x) args) args' [], 1)]]) (Plus [Term 1 []]))
+
+evalExpr env (SymbolicTensorExpr args sizeExpr name) = do
+  args' <- mapM (\arg -> evalExprDeep env arg) args
+  size' <- evalExpr env sizeExpr
+  size'' <- collectionToList size'
+  ns <- (mapM fromEgison size'') :: EgisonM [Integer]
+  let xs = map (\ms -> Value $ ScalarData (Div (Plus [Term 1 [(FunctionData (Just $ symbolScalarData "" (name ++ concat (map (\m -> "_" ++ m) (map show ms)))) (map (\x -> symbolScalarData "" $ show x) args) args' [], 1)]]) (Plus [Term 1 []])))
+               (map (\ms -> map toEgison ms) (enumTensorIndices ns))
+  fromTensor (Tensor ns (V.fromList xs) [])
+
 evalExpr env (IfExpr test expr expr') = do
   test <- evalExpr env test >>= fromWHNF
   evalExpr env $ if test then expr else expr'
@@ -384,7 +436,9 @@
  where
   extractBindings :: BindingExpr -> EgisonM [Binding]
   extractBindings ([name], expr) =
-    makeBindings [name] . (:[]) <$> newObjectRef env expr
+    case expr of
+      FunctionExpr args -> let Env frame _ = env in makeBindings [name] . (:[]) <$> newObjectRef (Env frame (Just $ varToVarWithIndices name)) expr
+      _ -> makeBindings [name] . (:[]) <$> newObjectRef env expr
   extractBindings (names, expr) =
     makeBindings names <$> (evalExpr env expr >>= fromTuple)
 
@@ -392,7 +446,7 @@
   let bindings' = evalState (concat <$> mapM extractBindings bindings) 0
   in recursiveBind env bindings' >>= flip evalExpr expr 
  where
-  extractBindings :: BindingExpr -> State Int [(String, EgisonExpr)]
+  extractBindings :: BindingExpr -> State Int [(Var, EgisonExpr)]
   extractBindings ([name], expr) = return [(name, expr)]
   extractBindings (names, expr) = do
     var <- genVar
@@ -401,12 +455,12 @@
         matcher = TupleExpr $ replicate k SomethingExpr
         nth n =
           let pattern = TuplePat $ flip map [1..k] $ \i ->
-                if i == n then PatVar "#_" else WildCard
-          in MatchExpr target matcher [(pattern, VarExpr "#_")]
+                if i == n then PatVar (stringToVar "#_") else WildCard
+          in MatchExpr target matcher [(pattern, VarExpr $ stringToVar "#_")]
     return ((var, expr) : map (second nth) (zip names [1..]))
 
-  genVar :: State Int String
-  genVar = modify (1+) >> gets (('#':) . show)
+  genVar :: State Int Var
+  genVar = modify (1+) >> gets (stringToVar . ('#':) . show)
 
 evalExpr env (TransposeExpr vars expr) = do
   syms <- evalExpr env vars >>= collectionToList
@@ -434,7 +488,7 @@
 evalExpr env (WithSymbolsExpr vars expr) = do
   symId <- fresh
   syms <- mapM (\var -> (newEvaluatedObjectRef (Value (symbolScalarData symId var)))) vars
-  let bindings = zip vars syms
+  let bindings = zip (map stringToVar vars) syms
   whnf <- evalExpr (extendEnv env bindings) expr
   case whnf of
     (Value (TensorData (Tensor ns xs js))) -> do
@@ -453,6 +507,9 @@
   isTmpSymbol symId (SupSubscript (ScalarData (Div (Plus [Term 1 [(Symbol id name is,n)]]) (Plus [Term 1 []]))))
     | symId == id = True
     | otherwise = False
+  isTmpSymbol symId (Userscript (ScalarData (Div (Plus [Term 1 [(Symbol id name is,n)]]) (Plus [Term 1 []]))))
+    | symId == id = True
+    | otherwise = False
   removeTmpscripts :: String -> WHNFData -> EgisonM WHNFData
   removeTmpscripts symId (Intermediate (ITensor (Tensor s xs is))) = do
     let (ds, js) = partition (isTmpSymbol symId) is
@@ -466,12 +523,12 @@
     
 
 evalExpr env (DoExpr bindings expr) = return $ Value $ IOFunc $ do
-  let body = foldr genLet (ApplyExpr expr $ TupleExpr [VarExpr "#1"]) bindings
+  let body = foldr genLet (ApplyExpr expr $ TupleExpr [VarExpr $ stringToVar "#1"]) bindings
   applyFunc env (Value $ Func Nothing env ["#1"] body) $ Value World
  where
   genLet (names, expr) expr' =
-    LetExpr [(["#1", "#2"], ApplyExpr expr $ TupleExpr [VarExpr "#1"])] $
-    LetExpr [(names, VarExpr "#2")] expr'
+    LetExpr [(map stringToVar ["#1", "#2"], ApplyExpr expr $ TupleExpr [VarExpr $ stringToVar "#1"])] $
+    LetExpr [(names, VarExpr $ stringToVar "#2")] expr'
 
 evalExpr env (IoExpr expr) = do
   io <- evalExpr env expr
@@ -533,16 +590,9 @@
     _ -> applyFunc env func (Value (makeTuple args))
 
 evalExpr env (ApplyExpr func arg) = do
-  func <- evalExpr env func >>= appendDFscripts 0
+  func <- evalExpr env func
   arg <- evalExpr env arg
---  arg <- evalExpr env arg >>= fromTupleWHNF
---  let k = fromIntegral (length arg)
---  arg <-  mapM (\(_,j) -> appendDFscripts 0 j) (zip [1..k] arg) >>= makeITuple
   case func of
-    Value (TensorData t@(Tensor ns fs js)) -> do
-      tMap (\f -> applyFunc env (Value f) arg >>= evalWHNF) t >>= fromTensor >>= return . Value >>= removeDFscripts
-    Intermediate (ITensor t@(Tensor ns fs js)) -> do
-      tMap (\f -> applyFunc env f arg) t >>= fromTensor
     Value (MemoizedFunc name ref hashRef env names body) -> do
       indices <- evalWHNF arg
       indices' <- mapM fromEgison $ fromTupleValue indices
@@ -558,6 +608,32 @@
           writeObjectRef ref (Value (MemoizedFunc name ref hashRef env names body))
           return whnf
     _ -> applyFunc env func arg >>= removeDFscripts
+-- evalExpr env (ApplyExpr func arg) = do
+--   func <- evalExpr env func >>= appendDFscripts 0
+--   arg <- evalExpr env arg
+-- --  arg <- evalExpr env arg >>= fromTupleWHNF
+-- --  let k = fromIntegral (length arg)
+-- --  arg <-  mapM (\(_,j) -> appendDFscripts 0 j) (zip [1..k] arg) >>= makeITuple
+--   case func of
+--     Value (TensorData t@(Tensor ns fs js)) -> do
+--       tMap (\f -> applyFunc env (Value f) arg >>= evalWHNF) t >>= fromTensor >>= return . Value >>= removeDFscripts
+--     Intermediate (ITensor t@(Tensor ns fs js)) -> do
+--       tMap (\f -> applyFunc env f arg) t >>= fromTensor
+--     Value (MemoizedFunc name ref hashRef env names body) -> do
+--       indices <- evalWHNF arg
+--       indices' <- mapM fromEgison $ fromTupleValue indices
+--       hash <- liftIO $ readIORef hashRef
+--       case HL.lookup indices' hash of
+--         Just objRef -> do
+--           evalRef objRef
+--         Nothing -> do
+--           whnf <- applyFunc env (Value (Func Nothing env names body)) arg
+--           retRef <- newEvaluatedObjectRef whnf
+--           hash <- liftIO $ readIORef hashRef
+--           liftIO $ writeIORef hashRef (HL.insert indices' retRef hash)
+--           writeObjectRef ref (Value (MemoizedFunc name ref hashRef env names body))
+--           return whnf
+--     _ -> applyFunc env func arg >>= removeDFscripts
 
 evalExpr env (WedgeApplyExpr func arg) = do
   func <- evalExpr env func >>= appendDFscripts 0
@@ -602,7 +678,7 @@
 
 evalExpr env (MatcherBFSExpr info) = return $ Value $ UserMatcher env BFSMode info
 evalExpr env (MatcherDFSExpr info) = return $ Value $ UserMatcher env DFSMode info
-
+ 
 evalExpr env (GenerateArrayExpr fnExpr (fstExpr, lstExpr)) = do
   fN <- (evalExpr env fstExpr >>= fromWHNF) :: EgisonM Integer
   eN <- (evalExpr env lstExpr >>= fromWHNF) :: EgisonM Integer
@@ -616,10 +692,20 @@
   size' <- evalExpr env sizeExpr
   size'' <- collectionToList size'
   ns <- (mapM fromEgison size'') :: EgisonM [Integer]
-  fn <- evalExpr env fnExpr
-  xs <-  mapM (\ms -> applyFunc env fn (Value (makeTuple ms))) (map (\ms -> map toEgison ms) (enumTensorIndices ns))
-  case (ns, xs) of
-    _ -> fromTensor (Tensor ns (V.fromList xs) [])
+  let Env frame maybe_vwi = env
+  xs <- mapM (\ms -> do
+    let env' = case maybe_vwi of
+                 Nothing -> env
+                 Just (VarWithIndices nameString indexList) -> Env frame (Just $ VarWithIndices nameString $ changeIndexList indexList ms)
+    fn <- evalExpr env' fnExpr
+    applyFunc env fn (Value (makeTuple ms))) 
+                (map (\ms -> map toEgison ms) (enumTensorIndices ns))
+  fromTensor (Tensor ns (V.fromList xs) [])
+ where 
+   changeIndexList :: [Index String] -> [EgisonValue] -> [Index String]
+   changeIndexList idxlist ms = map (\(i, m) -> case i of
+                                                  Superscript s -> Superscript (s ++ m)
+                                                  Subscript s -> Subscript (s ++ m)) $ zip idxlist (map show ms)
 
 evalExpr env (TensorContractExpr fnExpr tExpr) = do
   fn <- evalExpr env fnExpr
@@ -754,51 +840,74 @@
 evalWHNF (Intermediate (ITensor (Tensor ns whnfs js))) = do
   vals <- mapM evalWHNF (V.toList whnfs)
   return $ TensorData $ Tensor ns (V.fromList vals) js
---  vals <- mapM evalWHNF whnfs
---  return $ TensorData $ Tensor ns vals js
 evalWHNF coll = Collection <$> (fromCollection coll >>= fromMList >>= mapM evalRefDeep . Sq.fromList)
 
-fib :: Integer -> Integer
-fib 0 = 1
-fib 1 = 1
-fib n = (fib (n - 1)) + (fib (n - 2))
+addscript :: (Index EgisonValue, Tensor a) -> Tensor a
+addscript (subj, (Tensor s t i)) = (Tensor s t (i ++ [subj]))
 
+valuetoTensor2 :: WHNFData -> Tensor WHNFData
+valuetoTensor2 (Intermediate (ITensor t)) = t
+
 applyFunc :: Env -> WHNFData -> WHNFData -> EgisonM WHNFData
+applyFunc env (Value (TensorData (Tensor s1 t1 i1))) tds = do
+    tds <- fromTupleWHNF tds
+    if (length s1) > (length i1) && (all (\(Intermediate (ITensor (Tensor s u i))) -> ((length s) - (length i) == 1)) tds)
+       then do
+            symId <- fresh
+            let argnum = length tds
+                subjs = map (\symName -> Subscript $ symbolScalarData symId (show symName)) [1 .. argnum]
+                supjs = map (\symName -> Superscript $ symbolScalarData symId (show symName)) [1 .. argnum]
+            dot <- evalExpr env (VarExpr $ stringToVar ".")
+            makeITuple ((Value (TensorData (Tensor s1 t1 (i1 ++ supjs)))):(map Intermediate (map (ITensor . addscript) (zip subjs $ map valuetoTensor2 tds)))) >>= applyFunc env dot
+       else throwError $ Default "applyfunc"
+
+applyFunc env (Intermediate (ITensor (Tensor s1 t1 i1))) tds = do
+    tds <- fromTupleWHNF tds
+    if (length s1) > (length i1) && (all (\(Intermediate (ITensor (Tensor s u i))) -> ((length s) - (length i) == 1)) tds)
+       then do
+            symId <- fresh
+            let argnum = length tds
+                subjs = map (\symName -> Subscript $ symbolScalarData symId (show symName)) [1 .. argnum]
+                supjs = map (\symName -> Superscript $ symbolScalarData symId (show symName)) [1 .. argnum]
+            dot <- evalExpr env (VarExpr $ stringToVar ".")
+            makeITuple (map Intermediate (ITensor (Tensor s1 t1 (i1 ++ supjs)):(map (ITensor . addscript) (zip subjs $ map valuetoTensor2 tds)))) >>= applyFunc env dot
+       else throwError $ Default "applyfunc"
+
 applyFunc _ (Value (PartialFunc env n body)) arg = do
   refs <- fromTuple arg
   if n == fromIntegral (length refs)
-    then evalExpr (extendEnv env $ makeBindings (map (\n -> "::" ++ show n) [1..n]) refs) body
+    then evalExpr (extendEnv env $ makeBindings (map (\n -> (stringToVar $ "::" ++ show n)) [1..n]) refs) body
     else throwError $ ArgumentsNumWithNames ["partial"] (fromIntegral n) (length refs)
 applyFunc _ (Value (Func _ env [name] body)) arg = do
   ref <- newEvaluatedObjectRef arg
-  evalExpr (extendEnv env $ makeBindings [name] [ref]) body
+  evalExpr (extendEnv env $ makeBindings' [name] [ref]) body
 applyFunc _ (Value (Func _ env names body)) arg = do
   refs <- fromTuple arg
   if length names == length refs
-    then evalExpr (extendEnv env $ makeBindings names refs) body
+    then evalExpr (extendEnv env $ makeBindings' names refs) body
     else throwError $ ArgumentsNumWithNames names (length names) (length refs)
 applyFunc _ (Value (Proc _ env [name] body)) arg = do
   ref <- newEvaluatedObjectRef arg
-  evalExpr (extendEnv env $ makeBindings [name] [ref]) body
+  evalExpr (extendEnv env $ makeBindings' [name] [ref]) body
 applyFunc _ (Value (Proc _ env names body)) arg = do
   refs <- fromTuple arg
   if length names == length refs
-    then evalExpr (extendEnv env $ makeBindings names refs) body
+    then evalExpr (extendEnv env $ makeBindings' names refs) body
     else throwError $ ArgumentsNumWithNames names (length names) (length refs)
 applyFunc _ (Value (CFunc _ env name body)) arg = do
   refs <- fromTuple arg
   seqRef <- liftIO . newIORef $ Sq.fromList (map IElement refs)
   col <- liftIO . newIORef $ WHNF $ Intermediate $ ICollection $ seqRef
   if length refs > 0
-    then evalExpr (extendEnv env $ makeBindings [name] [col]) body
+    then evalExpr (extendEnv env $ makeBindings' [name] [col]) body
     else throwError $ ArgumentsNumWithNames [name] 1 0
 applyFunc env (Value (Macro [name] body)) arg = do
   ref <- newEvaluatedObjectRef arg
-  evalExpr (extendEnv env $ makeBindings [name] [ref]) body
+  evalExpr (extendEnv env $ makeBindings' [name] [ref]) body
 applyFunc env (Value (Macro names body)) arg = do
   refs <- fromTuple arg
   if length names == length refs
-    then evalExpr (extendEnv env $ makeBindings names refs) body
+    then evalExpr (extendEnv env $ makeBindings' names refs) body
     else throwError $ ArgumentsNumWithNames names (length names) (length refs)
 applyFunc _ (Value (PrimitiveFunc _ func)) arg = func arg
 applyFunc _ (Value (IOFunc m)) arg = do
@@ -813,10 +922,6 @@
   args <- tupleToList arg
   mExprs <- mapM extractScalar args
   return (Value (ScalarData (Div (Plus [(Term 1 [(Apply fn mExprs, 1)])]) (Plus [(Term 1 [])]))))
-applyFunc _ (Value fn@(UserIndexedData _ _)) arg = do
-  args <- tupleToList arg
-  mExprs <- mapM extractScalar args
-  return (Value (ScalarData (Div (Plus [(Term 1 [(Apply fn mExprs, 1)])]) (Plus [(Term 1 [])]))))
 applyFunc _ whnf _ = throwError $ TypeMismatch "function" whnf
 
 refArray :: WHNFData -> [EgisonValue] -> EgisonM WHNFData
@@ -902,15 +1007,18 @@
 newEvaluatedObjectRef :: WHNFData -> EgisonM ObjectRef
 newEvaluatedObjectRef = liftIO . newIORef . WHNF
 
-makeBindings :: [String] -> [ObjectRef] -> [Binding]
+makeBindings :: [Var] -> [ObjectRef] -> [Binding]
 makeBindings = zip
 
-recursiveBind :: Env -> [(String, EgisonExpr)] -> EgisonM Env
+makeBindings' :: [String] -> [ObjectRef] -> [Binding]
+makeBindings' xs = zip (map stringToVar xs)
+
+recursiveBind :: Env -> [(Var, EgisonExpr)] -> EgisonM Env
 recursiveBind env bindings = do
   let (names, exprs) = unzip bindings
   refs <- replicateM (length bindings) $ newObjectRef nullEnv UndefinedExpr
   let env' = extendEnv env $ makeBindings names refs
-  zipWithM_ (\ref (name,expr) ->
+  zipWithM_ (\ref (name,expr) -> do
                case expr of
                  MemoizedLambdaExpr names body -> do
                    hashRef <- liftIO $ newIORef HL.empty
@@ -923,14 +1031,20 @@
                    whnf <- evalExpr env' expr
                    case whnf of
                      (Value (CFunc _ env arg body)) -> liftIO . writeIORef ref . WHNF $ (Value (CFunc (Just name) env arg body))
+                 FunctionExpr args -> do
+                   let Env frame _ = env'
+                   liftIO . writeIORef ref . Thunk $ evalExpr (Env frame (Just $ varToVarWithIndices name)) $ FunctionExpr args
+                 GenerateTensorExpr _ _ -> do
+                   let Env frame _ = env'
+                   liftIO . writeIORef ref . Thunk $ evalExpr (Env frame (Just $ varToVarWithIndices name)) $ expr
                  _ -> liftIO . writeIORef ref . Thunk $ evalExpr env' expr)
             refs bindings
   return env'
 
-recursiveRebind :: Env -> (String, EgisonExpr) -> EgisonM Env
+recursiveRebind :: Env -> (Var, EgisonExpr) -> EgisonM Env
 recursiveRebind env (name, expr) = do
   case refVar env name of
-    Nothing -> throwError $ UnboundVariable name
+    Nothing -> throwError $ UnboundVariable $ show name
     Just ref -> case expr of
                   MemoizedLambdaExpr names body -> do
                     hashRef <- liftIO $ newIORef HL.empty
@@ -1227,7 +1341,7 @@
 primitivePatPatternMatch _ PPPatVar pattern = return ([pattern], [])
 primitivePatPatternMatch env (PPValuePat name) (ValuePat expr) = do
   ref <- lift $ newObjectRef env expr
-  return ([], [(name, ref)])
+  return ([], [(stringToVar name, ref)])
 primitivePatPatternMatch env (PPInductivePat name patterns) (InductivePat name' exprs)
   | name == name' && length patterns == length exprs =
     (concat *** concat) . unzip <$> zipWithM (primitivePatPatternMatch env) patterns exprs
@@ -1238,7 +1352,7 @@
 primitiveDataPatternMatch PDWildCard _ = return []
 primitiveDataPatternMatch (PDPatVar name) whnf = do
   ref <- lift $ newEvaluatedObjectRef whnf
-  return [(name, ref)]
+  return [(stringToVar name, ref)]
 primitiveDataPatternMatch (PDInductivePat name patterns) whnf = do
   case whnf of
     Intermediate (IInductiveData name' refs) | name == name' -> do
diff --git a/hs-src/Language/Egison/Desugar.hs b/hs-src/Language/Egison/Desugar.hs
--- a/hs-src/Language/Egison/Desugar.hs
+++ b/hs-src/Language/Egison/Desugar.hs
@@ -57,7 +57,7 @@
 
 desugar :: EgisonExpr -> DesugarM EgisonExpr
 desugar (AlgebraicDataMatcherExpr patterns) = do
-  matcherName <- fresh
+  matcherName <- freshV
   matcherRef <- return $ VarExpr matcherName
   matcher <- genMatcherClauses patterns matcherRef
   return $ LetRecExpr [([matcherName], matcher)] matcherRef
@@ -74,7 +74,7 @@
       genMainClause patterns matcher = do
         clauses <- genClauses patterns
         return (PPValuePat "val", TupleExpr []
-               ,[(PDPatVar "tgt", (MatchExpr (TupleExpr [(VarExpr "val"), (VarExpr "tgt")]) 
+               ,[(PDPatVar "tgt", (MatchExpr (TupleExpr [(VarExpr $ stringToVar "val"), (VarExpr $ stringToVar "tgt")]) 
                                              (TupleExpr [matcher, matcher]) 
                                              clauses))])
         where
@@ -89,7 +89,7 @@
         
           genMatchingPattern :: (String, [EgisonExpr]) -> DesugarM (EgisonPattern, EgisonPattern)
           genMatchingPattern (name, patterns) = do
-            names <- mapM (const fresh) patterns
+            names <- mapM (const freshV) patterns
             return $ ((InductivePat name (map PatVar names))  
                      ,(InductivePat name (map (ValuePat . VarExpr) names)))
           
@@ -107,8 +107,8 @@
           
           genPrimitiveDataPat :: (String, [EgisonExpr]) -> DesugarM (PrimitiveDataPattern, [EgisonExpr])
           genPrimitiveDataPat (name, patterns) = do
-            patterns' <- mapM (const fresh) patterns 
-            return (PDInductivePat (capitalize name) $ map PDPatVar patterns', map VarExpr patterns')
+            patterns' <- mapM (const freshV) patterns 
+            return (PDInductivePat (capitalize name) $ map PDPatVar $ map show patterns', map VarExpr patterns')
 
           capitalize :: String -> String
           capitalize (x:xs) = toUpper x : xs
@@ -116,7 +116,7 @@
       
       genSomethingClause :: DesugarM (PrimitivePatPattern, EgisonExpr, [(PrimitiveDataPattern, EgisonExpr)])
       genSomethingClause = 
-        return (PPPatVar, (TupleExpr [SomethingExpr]), [(PDPatVar "tgt", CollectionExpr [ElementExpr (VarExpr "tgt")])])
+        return (PPPatVar, (TupleExpr [SomethingExpr]), [(PDPatVar "tgt", CollectionExpr [ElementExpr (VarExpr $ stringToVar "tgt")])])
     
       matchingSuccess :: EgisonExpr
       matchingSuccess = CollectionExpr [ElementExpr $ TupleExpr []]
@@ -126,30 +126,55 @@
 
 desugar (MatchAllLambdaExpr matcher clause) = do
   name <- fresh
-  desugar $ LambdaExpr [TensorArg name] (MatchAllExpr (VarExpr name) matcher clause)
+  desugar $ LambdaExpr [TensorArg name] (MatchAllExpr (VarExpr $ stringToVar name) matcher clause)
 
 desugar (MatchLambdaExpr matcher clauses) = do
   name <- fresh
-  desugar $ LambdaExpr [TensorArg name] (MatchExpr (VarExpr name) matcher clauses)
+  desugar $ LambdaExpr [TensorArg name] (MatchExpr (VarExpr $ stringToVar name) matcher clauses)
 
 desugar (ArrayRefExpr expr nums) =
   case nums of
-    (TupleExpr nums') -> desugar $ IndexedExpr False expr (map Subscript nums')
-    _ -> desugar $ IndexedExpr False expr [Subscript nums]
-  
-desugar (IndexedExpr b expr indices) = 
-  IndexedExpr b <$> desugar expr <*> (mapM desugarIndex indices)
+    (TupleExpr nums') -> desugar $ IndexedExpr True expr (map Subscript nums')
+    _ -> desugar $ IndexedExpr True expr [Subscript nums]
 
-desugar (SubrefsExpr expr1 expr2) = 
-  SubrefsExpr <$> desugar expr1 <*> desugar expr2
+desugar (IndexedExpr b expr indices)
+  | endWithThreeDots expr = case expr of
+                              (VarExpr name) -> let x = show name in desugar $ IndexedExpr False (VarExpr $ stringToVar $ take ((length x)-3) x) indices
+  | otherwise = case indices of
+                 [MultiSubscript x y] -> case (x, y) of
+                                           (IntegerExpr _, IntegerExpr _) -> return $ SubrefsExpr b expr (ApplyExpr (VarExpr $ stringToVar "between") (TupleExpr [x, y]))
+                                           (TupleExpr [IndexedExpr b1 e1 [n1]], TupleExpr [IndexedExpr b2 e2 [n2]]) -> do
+                                             k <- fresh
+                                             return $ SubrefsExpr b expr (ApplyExpr (VarExpr $ stringToVar "map")
+                                                                                    (TupleExpr [(LambdaExpr [TensorArg k] (IndexedExpr b1 e1 [(Subscript $ VarExpr $ stringToVar k)])),
+                                                                                                (ApplyExpr (VarExpr $ stringToVar "between") (TupleExpr [(fromIndexToExpr n1), (fromIndexToExpr n2)]))]))
+                 [MultiSuperscript x y] -> case (x, y) of
+                                             (IntegerExpr _, IntegerExpr _) -> return $ SubrefsExpr b expr (ApplyExpr (VarExpr $ stringToVar "between") (TupleExpr [x, y]))
+                                             (TupleExpr [IndexedExpr b1 e1 [n1]], TupleExpr [IndexedExpr b2 e2 [n2]]) -> do    
+                                               k <- fresh
+                                               return $ SuprefsExpr b expr (ApplyExpr (VarExpr $ stringToVar "map")
+                                                                                      (TupleExpr [(LambdaExpr [TensorArg k] (IndexedExpr b1 e1 [(Subscript $ VarExpr $ stringToVar k)])),
+                                                                                                  (ApplyExpr (VarExpr $ stringToVar "between") (TupleExpr [(fromIndexToExpr n1), (fromIndexToExpr n2)]))]))
+                 _ -> IndexedExpr b <$> desugar expr <*> (mapM desugarIndex indices)
+ where
+  endWithThreeDots :: EgisonExpr -> Bool
+  endWithThreeDots (VarExpr name) = (take 3 $ reverse (show name)) == "..."
+  endWithThreeDots _ = False
+  fromIndexToExpr :: Index EgisonExpr -> EgisonExpr
+  fromIndexToExpr (Subscript a) = a
+  fromIndexToExpr (Superscript a) = a
+  fromIndexToExpr (SupSubscript a) = a
 
-desugar (SuprefsExpr expr1 expr2) = 
-  SuprefsExpr <$> desugar expr1 <*> desugar expr2
+desugar (SubrefsExpr bool expr1 expr2) = 
+  SubrefsExpr bool <$> desugar expr1 <*> desugar expr2
 
+desugar (SuprefsExpr bool expr1 expr2) = 
+  SuprefsExpr bool <$> desugar expr1 <*> desugar expr2
+
 desugar (PowerExpr expr1 expr2) = do
   expr1' <- desugar expr1
   expr2' <- desugar expr2
-  return $ ApplyExpr (VarExpr "**") (TupleExpr [expr1', expr2'])
+  return $ ApplyExpr (VarExpr $ stringToVar "**") (TupleExpr [expr1', expr2'])
 
 desugar (ArrayBoundsExpr expr) = do
   expr' <- desugar expr
@@ -187,24 +212,37 @@
 desugar (LambdaExpr names expr) = do
   let (rtnames, rhnames) = span (\name -> case name of
                                             TensorArg _ -> True
-                                            InvertedScalarArg _ -> False
-                                            ScalarArg _ -> False) (reverse names)
+                                            _ -> False) (reverse names)
   case rhnames of
     [] -> do expr' <- desugar expr
              return $ LambdaExpr names expr'
+
     (InvertedScalarArg rhname:rhnames') -> do
-      desugar $ LambdaExpr (reverse rhnames' ++ [TensorArg rhname] ++ reverse rtnames)
-                  (TensorMapExpr (LambdaExpr [TensorArg rhname] expr) (FlipIndicesExpr (VarExpr rhname)))
+      let (rtnames2, rhnames2) = span (\name -> case name of
+                                                _ -> False) rhnames'
+      case rhnames2 of
+        [] -> desugar $ LambdaExpr (reverse rhnames' ++ [TensorArg rhname] ++ reverse rtnames)
+                          (TensorMapExpr (LambdaExpr [TensorArg rhname] expr) (FlipIndicesExpr (VarExpr $ stringToVar rhname)))
+        (ScalarArg rhname2:rhnames2') ->
+          desugar $ LambdaExpr (reverse rhnames2' ++ [TensorArg rhname2] ++ rtnames2 ++ [TensorArg rhname] ++ reverse rtnames)
+                      (TensorMap2Expr (LambdaExpr [TensorArg rhname2, TensorArg rhname] expr) (VarExpr $ stringToVar rhname2) (FlipIndicesExpr (VarExpr $ stringToVar rhname)))
+        (InvertedScalarArg rhname2:rhnames2') ->
+          desugar $ LambdaExpr (reverse rhnames2' ++ [TensorArg rhname2] ++ rtnames2 ++ [TensorArg rhname] ++ reverse rtnames)
+                      (TensorMap2Expr (LambdaExpr [TensorArg rhname2, TensorArg rhname] expr) (FlipIndicesExpr (VarExpr $ stringToVar rhname2)) (FlipIndicesExpr (VarExpr $ stringToVar rhname)))
+
     (ScalarArg rhname:rhnames') -> do
       let (rtnames2, rhnames2) = span (\name -> case name of
                                                 TensorArg _ -> True
-                                                ScalarArg _ -> False) rhnames'
+                                                _ -> False) rhnames'
       case rhnames2 of
         [] -> desugar $ LambdaExpr (reverse rhnames' ++ [TensorArg rhname] ++ reverse rtnames)
-                          (TensorMapExpr (LambdaExpr [TensorArg rhname] expr) (VarExpr rhname))
+                          (TensorMapExpr (LambdaExpr [TensorArg rhname] expr) (VarExpr $ stringToVar rhname))
         (ScalarArg rhname2:rhnames2') ->
           desugar $ LambdaExpr (reverse rhnames2' ++ [TensorArg rhname2] ++ rtnames2 ++ [TensorArg rhname] ++ reverse rtnames)
-                      (TensorMap2Expr (LambdaExpr [TensorArg rhname2, TensorArg rhname] expr) (VarExpr rhname2) (VarExpr rhname))
+                      (TensorMap2Expr (LambdaExpr [TensorArg rhname2, TensorArg rhname] expr) (VarExpr $ stringToVar rhname2) (VarExpr $ stringToVar rhname))
+        (InvertedScalarArg rhname2:rhnames2') ->
+          desugar $ LambdaExpr (reverse rhnames2' ++ [TensorArg rhname2] ++ rtnames2 ++ [TensorArg rhname] ++ reverse rtnames)
+                      (TensorMap2Expr (LambdaExpr [TensorArg rhname2, TensorArg rhname] expr) (FlipIndicesExpr (VarExpr $ stringToVar rhname2)) (VarExpr $ stringToVar rhname))
 
 desugar (MemoizedLambdaExpr names expr) = do
   expr' <- desugar expr
@@ -348,11 +386,11 @@
   return $ CApplyExpr expr0' expr1'
 
 desugar (VarExpr name) = do
-  asks $ maybe (VarExpr name) id . lookup name
+  asks $ maybe (VarExpr name) id . lookup (show name)
 
 desugar FreshVarExpr = do
   id <- fresh
-  return (VarExpr (":::" ++ id))
+  return (VarExpr $ stringToVar (":::" ++ id))
 
 desugar (MatcherBFSExpr matcherInfo) = do
   matcherInfo' <- desugarMatcherInfo matcherInfo
@@ -366,15 +404,15 @@
 
 desugar (PartialExpr n expr) = do
   expr' <- desugar expr
-  return $ LetRecExpr [(["::0"], PartialExpr n expr')] (VarExpr "::0")
+  return $ LetRecExpr [([stringToVar "::0"], PartialExpr n expr')] (VarExpr $ stringToVar "::0")
 
 desugar (QuoteExpr expr) = do
   expr' <- desugar expr
   return $ QuoteExpr expr'
 
-desugar (QuoteFunctionExpr expr) = do
+desugar (QuoteSymbolExpr expr) = do
   expr' <- desugar expr
-  return $ QuoteFunctionExpr expr'
+  return $ QuoteSymbolExpr expr'
 
 desugar (WedgeExpr (ApplyExpr expr0 expr1)) = do
   expr0' <- desugar expr0
@@ -387,6 +425,7 @@
 desugarIndex (Superscript expr) = desugar expr >>= return . Superscript
 desugarIndex (Subscript expr) = desugar expr >>= return . Subscript
 desugarIndex (SupSubscript expr) = desugar expr >>= return . SupSubscript
+desugarIndex (Userscript expr) = desugar expr >>= return . Userscript
 
 desugarPattern :: EgisonPattern -> DesugarM EgisonPattern
 desugarPattern pattern = LetPat (map makeBinding $ S.elems $ collectName pattern) <$> desugarPattern' pattern 
@@ -402,9 +441,8 @@
    collectName (PApplyPat _ patterns) = collectNames patterns
    collectName (DApplyPat _ patterns) = collectNames patterns
    collectName (LoopPat _ (LoopRange _ _ endNumPat) pattern1 pattern2) = collectName endNumPat `S.union` collectName pattern1 `S.union` collectName pattern2
---   collectName (LoopPat _ (LoopRange _ _ endNumPat) pattern1 pattern2) = collectName pattern1 `S.union` collectName pattern2
    collectName (LetPat _ pattern) = collectName pattern
-   collectName (IndexedPat (PatVar name) _) = S.singleton name
+   collectName (IndexedPat (PatVar name) _) = S.singleton $ show name
    collectName (OrPat patterns) = collectNames patterns
    collectName (DivPat pattern1 pattern2) = collectName pattern1 `S.union` collectName pattern2
    collectName (PlusPat patterns) = collectNames patterns
@@ -413,7 +451,7 @@
    collectName _ = S.empty
    
    makeBinding :: String -> BindingExpr
-   makeBinding name = ([name], HashExpr [])
+   makeBinding name = ([stringToVar name], HashExpr [])
 
 desugarPattern' :: EgisonPattern -> DesugarM EgisonPattern
 desugarPattern' (ValuePat expr) = ValuePat <$> desugar expr
diff --git a/hs-src/Language/Egison/MathOutput.hs b/hs-src/Language/Egison/MathOutput.hs
--- a/hs-src/Language/Egison/MathOutput.hs
+++ b/hs-src/Language/Egison/MathOutput.hs
@@ -27,7 +27,7 @@
                           Left err -> input
                           Right val -> "#latex|" ++ showMathExprLatex val ++ "|#"
 
-data MathExpr = Atom String
+data MathExpr = Atom String [MathIndex]
               | NegativeAtom String
               | Plus [MathExpr]
               | Multiply [MathExpr]
@@ -38,7 +38,7 @@
               | Collection [MathExpr]
               | Exp MathExpr
               | Quote MathExpr
-              | Partial String [String]
+              | Partial MathExpr [MathExpr]
               deriving (Eq, Show)
 
 data MathIndex = Super MathExpr
@@ -54,7 +54,7 @@
 showMathIndexAsciiMath (Sub a) = showMathExprAsciiMath a
 
 showMathExprAsciiMath :: MathExpr -> String
-showMathExprAsciiMath (Atom func) = func
+showMathExprAsciiMath (Atom func []) = func
 showMathExprAsciiMath (NegativeAtom func) = "-" ++ func
 showMathExprAsciiMath (Plus []) = ""
 showMathExprAsciiMath (Plus (x:xs)) = showMathExprAsciiMath x ++ showMathExprAsciiMathForPlus xs
@@ -62,7 +62,7 @@
   showMathExprAsciiMathForPlus :: [MathExpr] -> String
   showMathExprAsciiMathForPlus [] = ""
   showMathExprAsciiMathForPlus ((NegativeAtom a):xs) = " - " ++ a ++ showMathExprAsciiMathForPlus xs
-  showMathExprAsciiMathForPlus ((Multiply (NegativeAtom a:ys)):xs) = " - " ++ showMathExprAsciiMath (Multiply ((Atom a):ys)) ++ " " ++ showMathExprAsciiMathForPlus xs
+  showMathExprAsciiMathForPlus ((Multiply (NegativeAtom a:ys)):xs) = " - " ++ showMathExprAsciiMath (Multiply ((Atom a []):ys)) ++ " " ++ showMathExprAsciiMathForPlus xs
   showMathExprAsciiMathForPlus (x:xs) = showMathExprAsciiMath x ++ " + " ++ showMathExprAsciiMathForPlus xs
 showMathExprAsciiMath (Multiply []) = ""
 showMathExprAsciiMath (Multiply [a]) = showMathExprAsciiMath a
@@ -70,7 +70,7 @@
 showMathExprAsciiMath (Multiply lvs) = showMathExprAsciiMath' (head lvs) ++ " " ++ showMathExprAsciiMath (Multiply (tail lvs))
 showMathExprAsciiMath (Power lv1 lv2) = showMathExprAsciiMath lv1 ++ "^" ++ showMathExprAsciiMath lv2
 showMathExprAsciiMath (Func f lvs) = case f of
-                                       Atom "/" -> if length lvs == 2 then "frac{" ++ showMathExprAsciiMath (head lvs) ++ "}{" ++ showMathExprAsciiMath (lvs !! 1) ++ "}"
+                                       Atom "/" [] -> if length lvs == 2 then "frac{" ++ showMathExprAsciiMath (head lvs) ++ "}{" ++ showMathExprAsciiMath (lvs !! 1) ++ "}"
                                                                         else showMathExprAsciiMath f ++ "(" ++ showMathExprAsciiMathArg lvs ++ ")"
                                        _ -> showMathExprAsciiMath f ++ "(" ++ showMathExprAsciiMathArg lvs ++ ")"
 showMathExprAsciiMath (Tensor lvs mis)
@@ -105,8 +105,19 @@
 --
 
 showMathExprLatex :: MathExpr -> String
-showMathExprLatex (Atom a) = a
-showMathExprLatex (Partial a is) = a ++ "_{" ++ concat is ++ "}"
+showMathExprLatex (Atom a []) = a
+showMathExprLatex (Atom a xs) = a ++ showMathExprLatexScript xs
+showMathExprLatex (Partial f xs) = "\\frac{" ++ convertToPartial (f, length xs) ++ "}{" ++ showPartial xs ++ "}"
+                                         where showPartial :: [MathExpr] -> String
+                                               showPartial xs = let lx = elemCount xs in convertToPartial2 (head lx) ++ foldr (\x acc -> " " ++ convertToPartial2 x ++ acc) "" (tail lx)
+
+                                               convertToPartial :: (MathExpr, Int) -> String
+                                               convertToPartial (x, 1) = "\\partial " ++ showMathExprLatex x
+                                               convertToPartial (x, n) = "\\partial^" ++ show n ++ " " ++ showMathExprLatex x
+
+                                               convertToPartial2 :: (MathExpr, Int) -> String
+                                               convertToPartial2 (x, 1) = "\\partial " ++ showMathExprLatex x
+                                               convertToPartial2 (x, n) = "\\partial " ++ showMathExprLatex x ++ "^"  ++ show n
 showMathExprLatex (NegativeAtom a) = "-" ++ a
 showMathExprLatex (Plus []) = ""
 showMathExprLatex (Plus (x:xs)) = showMathExprLatex x ++ showMathExprLatexForPlus xs
@@ -114,17 +125,17 @@
   showMathExprLatexForPlus :: [MathExpr] -> String
   showMathExprLatexForPlus [] = ""
   showMathExprLatexForPlus ((NegativeAtom a):xs) = " - " ++ a ++ showMathExprLatexForPlus xs
-  showMathExprLatexForPlus ((Multiply (NegativeAtom a:ys)):xs) = " - " ++ showMathExprLatex (Multiply ((Atom a):ys)) ++ showMathExprLatexForPlus xs
+  showMathExprLatexForPlus ((Multiply (NegativeAtom a:ys)):xs) = " - " ++ showMathExprLatex (Multiply ((Atom a []):ys)) ++ showMathExprLatexForPlus xs
   showMathExprLatexForPlus (x:xs) = " + " ++  showMathExprLatex x ++ showMathExprLatexForPlus xs
 showMathExprLatex (Multiply []) = ""
 showMathExprLatex (Multiply [x]) = showMathExprLatex x
-showMathExprLatex (Multiply (Atom "1":xs)) = showMathExprLatex (Multiply xs)
+showMathExprLatex (Multiply (Atom "1" []:xs)) = showMathExprLatex (Multiply xs)
 showMathExprLatex (Multiply (NegativeAtom "1":xs)) = "-" ++ showMathExprLatex (Multiply xs)
 showMathExprLatex (Multiply (x:xs)) = showMathExprLatex' x ++ " " ++ showMathExprLatex (Multiply xs)
 showMathExprLatex (Power lv1 lv2) = showMathExprLatex lv1 ++ "^" ++ showMathExprLatex lv2
-showMathExprLatex (Func (Atom "sqrt") [x]) = "\\sqrt{" ++ showMathExprLatex x ++ "}"
-showMathExprLatex (Func (Atom "rt") [x, y]) = "\\sqrt[" ++ showMathExprLatex x ++ "]{" ++ showMathExprLatex y ++ "}"
-showMathExprLatex (Func (Atom "/") [x, y]) = "\\frac{" ++ showMathExprLatex x ++ "}{" ++ showMathExprLatex y ++ "}"
+showMathExprLatex (Func (Atom "sqrt" []) [x]) = "\\sqrt{" ++ showMathExprLatex x ++ "}"
+showMathExprLatex (Func (Atom "rt" []) [x, y]) = "\\sqrt[" ++ showMathExprLatex x ++ "]{" ++ showMathExprLatex y ++ "}"
+showMathExprLatex (Func (Atom "/" []) [x, y]) = "\\frac{" ++ showMathExprLatex x ++ "}{" ++ showMathExprLatex y ++ "}"
 showMathExprLatex (Func f xs) = showMathExprLatex f ++ "(" ++ showMathExprLatexArg xs ", " ++ ")"
 showMathExprLatex (Tensor xs mis) = case head xs of
                                        Tensor _ _ -> "\\begin{pmatrix} " ++ showMathExprLatexVectors xs ++ "\\end{pmatrix}" ++ showMathExprLatexScript mis
@@ -144,12 +155,12 @@
 showMathExprLatexArg lvs s = showMathExprLatex (head lvs) ++ s ++ showMathExprLatexArg  (tail lvs) s
 
 showMathExprLatexSuper :: MathIndex -> String
-showMathExprLatexSuper (Super (Atom "#")) = "\\#"
+showMathExprLatexSuper (Super (Atom "#" [])) = "\\#"
 showMathExprLatexSuper (Super x) = showMathExprLatex x
 showMathExprLatexSuper (Sub x) = "\\;"
 
 showMathExprLatexSub :: MathIndex -> String
-showMathExprLatexSub (Sub (Atom "#")) = "\\#"
+showMathExprLatexSub (Sub (Atom "#" [])) = "\\#"
 showMathExprLatexSub (Sub x) = showMathExprLatex x
 showMathExprLatexSub (Super x) = "\\;"
 
@@ -180,8 +191,21 @@
     first <- letter <|> symbol <|> digit
     rest <- many (letter <|> digit <|> symbol)
     let atom = first : rest
-    option (Atom atom) $ do is <- many1 (char '|' >> many digit)
-                            return $ Partial atom is
+    ys <- many parseScript
+    return $ Atom atom ys
+
+parseAtom' :: Parser MathExpr
+parseAtom' = do 
+    first <- letter <|> symbol <|> digit
+    rest <- many (letter <|> digit <|> symbol)
+    let atom = first : rest
+    return $ Atom atom []
+
+parsePartial :: Parser MathExpr
+parsePartial = do
+    xs <- parseAtom
+    is <- many1 (char '|' >> parseAtom)
+    return $ Partial xs is
   
 parseNegativeAtom :: Parser MathExpr
 parseNegativeAtom = do
@@ -195,7 +219,7 @@
 parseList = sepEndBy parseExpr spaces
 
 parseScript :: Parser MathIndex
-parseScript = (Sub <$> (char '_' >> parseExpr)) <|> (Super <$> (char '~' >> parseExpr))
+parseScript = (Sub <$> (char '_' >> parseAtom')) <|> (Super <$> (char '~' >> parseAtom'))
 
 parsePlus :: Parser MathExpr
 parsePlus = do
@@ -262,6 +286,7 @@
 
 parseExpr' :: Parser MathExpr
 parseExpr' = parseNegativeAtom
+         <|> try parsePartial
          <|> parseAtom
          <|> parseQuote
          <|> try parseExp
@@ -276,3 +301,7 @@
 parseExpr = do
     x <- parseExpr'
     option x $ Power x <$> try (char '^' >> parseExpr')
+
+elemCount :: Eq a => [a] -> [(a, Int)]
+elemCount [] = []
+elemCount (x:xs) = (x, (length $ filter (== x) xs) + 1) : elemCount (filter (/= x) xs)
diff --git a/hs-src/Language/Egison/Parser.hs b/hs-src/Language/Egison/Parser.hs
--- a/hs-src/Language/Egison/Parser.hs
+++ b/hs-src/Language/Egison/Parser.hs
@@ -38,6 +38,7 @@
 import qualified Data.Set as Set
 import Data.Traversable (mapM)
 import Data.Ratio
+import Data.List.Split (splitOn)
 
 import Text.Parsec
 import Text.Parsec.String
@@ -135,18 +136,11 @@
       <?> "top-level expression"
 
 defineExpr :: Parser EgisonTopExpr
-defineExpr = try (parens (keywordDefine >> Define <$> varNameWithIndexType <*> expr))
+defineExpr = try (parens (keywordDefine >> Define <$> (char '$' >> identVar) <*> expr))
          <|> try (parens (do keywordDefine
-                             (VarWithIndices name is) <- varNameWithIndices
+                             (VarWithIndices name is) <- (char '$' >> identVarWithIndices)
                              body <- expr
                              return $ Define (Var name (map f is)) (WithSymbolsExpr (map g is) (TransposeExpr (CollectionExpr (map (ElementExpr . h) is)) body))))
---defineExpr = try 
---                 (do keywordDefine
---                     (VarWithIndices name is) <- varNameWithIndices
---                     body <- expr
---                     return $ Define (Var name (map f is)) (WithSymbolsExpr (map g is) (TransposeExpr (CollectionExpr (map (ElementExpr . h) is)) body)))
---            <|> (keywordDefine >> Define <$> varNameWithIndexType <*> expr)
---defineExpr = (keywordDefine >> Define <$> varNameWithIndexType <*> expr)
  where
   f (Superscript _) = Superscript ()
   f (Subscript _) = Subscript ()
@@ -154,12 +148,12 @@
   g (Superscript i) = i
   g (Subscript i) = i
   g (SupSubscript i) = i
-  h (Superscript i) = (VarExpr i)
-  h (Subscript i) = (VarExpr i)
-  h (SupSubscript i) = (VarExpr i)
+  h (Superscript i) = (VarExpr $ stringToVar i)
+  h (Subscript i) = (VarExpr $ stringToVar i)
+  h (SupSubscript i) = (VarExpr $ stringToVar i)
 
 redefineExpr :: Parser EgisonTopExpr
-redefineExpr = (keywordRedefine <|> keywordSet) >> Redefine <$> varNameWithIndexType <*> expr
+redefineExpr = (keywordRedefine <|> keywordSet) >> Redefine <$> (char '$' >> identVar) <*> expr
 
 testExpr :: Parser EgisonTopExpr
 testExpr = keywordTest >> Test <$> expr
@@ -178,12 +172,28 @@
 
 expr :: Parser EgisonExpr
 expr = P.lexeme lexer (do expr0 <- expr' <|> quoteExpr'
-                          expr1 <- option expr0 $ IndexedExpr False expr0 <$> many1 (try (char '_' >> expr' >>= return . Subscript)
-                                                                           <|> try (char '~' >> expr' >>= return . Superscript)
-                                                                           <|> try (string "~_" >> expr' >>= return . SupSubscript))
-                          expr2 <- option expr1 $ UserIndexedExpr expr1 <$> many1 (try $ char '|' >> expr' >>= return . Userscript)
-                          option expr2 $ PowerExpr expr1 <$> (try $ char '^' >> expr'))
-                          
+                          expr1 <- option expr0 $ try (string "..." >> IndexedExpr False expr0 <$> parseindex)
+                                                  <|> IndexedExpr True expr0 <$> parseindex
+                          option expr1 $ PowerExpr expr1 <$> (try $ char '^' >> expr'))
+                            where parseindex :: Parser [Index EgisonExpr]
+                                  parseindex = many1 (try (do
+                                                           char '_' 
+                                                           e1 <- expr'
+                                                           string "..._"
+                                                           e2 <- expr'
+                                                           return $ MultiSubscript e1 e2)
+                                                 <|> try (do
+                                                           char '~'
+                                                           e1 <- expr'
+                                                           string "...~"
+                                                           e2 <- expr'
+                                                           return $ MultiSuperscript e1 e2)
+                                                 <|> try (char '_' >> expr' >>= return . Subscript)
+                                                 <|> try (char '~' >> expr' >>= return . Superscript)
+                                                 <|> try (string "~_" >> expr' >>= return . SupSubscript)
+                                                 <|> try (char '|' >> expr' >>= return . Userscript))
+
+
 quoteExpr' :: Parser EgisonExpr
 quoteExpr' = char '\'' >> QuoteExpr <$> expr'
 
@@ -200,7 +210,7 @@
              <|> try hashExpr
              <|> collectionExpr
 --             <|> quoteExpr
-             <|> quoteFunctionExpr
+             <|> quoteSymbolExpr
              <|> wedgeExpr
              <|> parens (ifExpr
                          <|> lambdaExpr
@@ -235,6 +245,7 @@
                          <|> arrayBoundsExpr
                          <|> arrayRefExpr
                          <|> generateTensorExpr
+                         <|> symbolicTensorExpr
                          <|> tensorExpr
                          <|> tensorContractExpr
                          <|> tensorMapExpr
@@ -245,11 +256,13 @@
                          <|> pmapExpr
                          <|> subrefsExpr
                          <|> suprefsExpr
+                         <|> userrefsExpr
+                         <|> functionWithArgExpr
                          )
              <?> "expression")
 
 varExpr :: Parser EgisonExpr
-varExpr = VarExpr <$> ident
+varExpr = VarExpr <$> identVarWithoutIndex
 
 freshVarExpr :: Parser EgisonExpr
 freshVarExpr = char '#' >> return FreshVarExpr
@@ -293,9 +306,18 @@
 wedgeExpr :: Parser EgisonExpr
 wedgeExpr = char '!' >> WedgeExpr <$> expr
 
-quoteFunctionExpr :: Parser EgisonExpr
-quoteFunctionExpr = char '`' >> QuoteFunctionExpr <$> expr
+functionWithArgExpr :: Parser EgisonExpr
+functionWithArgExpr = keywordFunction >> FunctionExpr <$> (between lp rp $ sepEndBy expr whiteSpace)
+  where
+    lp = P.lexeme lexer (char '[')
+    rp = char ']'
 
+symbolicTensorExpr :: Parser EgisonExpr
+symbolicTensorExpr = keywordSymbolicTensor >> SymbolicTensorExpr <$> (brackets $ sepEndBy expr whiteSpace) <*> expr <*> ident
+
+quoteSymbolExpr :: Parser EgisonExpr
+quoteSymbolExpr = char '`' >> QuoteSymbolExpr <$> expr
+
 matchAllExpr :: Parser EgisonExpr
 matchAllExpr = keywordMatchAll >> MatchAllExpr <$> expr <*> expr <*> matchClause
 
@@ -399,7 +421,7 @@
 memoizeBinding = brackets $ (,,) <$> expr <*> expr <*> expr
 
 cambdaExpr :: Parser EgisonExpr
-cambdaExpr = keywordCambda >> CambdaExpr <$> varName <*> expr
+cambdaExpr = keywordCambda >> char '$' >> CambdaExpr <$> ident <*> expr
 
 procedureExpr :: Parser EgisonExpr
 procedureExpr = keywordProcedure >> ProcedureExpr <$> varNames <*> expr
@@ -423,7 +445,7 @@
 withSymbolsExpr = keywordWithSymbols >> WithSymbolsExpr <$> (braces $ sepEndBy ident whiteSpace) <*> expr
 
 doExpr :: Parser EgisonExpr
-doExpr = keywordDo >> DoExpr <$> statements <*> option (ApplyExpr (VarExpr "return") (TupleExpr [])) expr
+doExpr = keywordDo >> DoExpr <$> statements <*> option (ApplyExpr (VarExpr $ stringToVar "return") (TupleExpr [])) expr
 
 statements :: Parser [BindingExpr]
 statements = braces $ sepEndBy statement whiteSpace
@@ -437,36 +459,15 @@
 bindings = braces $ sepEndBy binding whiteSpace
 
 binding :: Parser BindingExpr
-binding = brackets $ (,) <$> varNames <*> expr
+binding = brackets $ (,) <$> varNames' <*> expr
 
 varNames :: Parser [String]
-varNames = return <$> varName
-            <|> brackets (sepEndBy varName whiteSpace) 
-
-varName :: Parser String
-varName = char '$' >> ident
-
-varNameWithIndexType :: Parser Var
-varNameWithIndexType = P.lexeme lexer (do
-  char '$'
-  name <- ident
-  is <- many indexType
-  return $ Var name is)
-
-indexType :: Parser (Index ())
-indexType = try (char '~' >> return (Superscript ()))
-        <|> try (char '_' >> return (Subscript ()))
-
-varNameWithIndices :: Parser VarWithIndices
-varNameWithIndices = P.lexeme lexer (do
-  char '$'
-  name <- ident
-  is <- many indexForVar
-  return $ VarWithIndices name is)
+varNames = return <$> (char '$' >> ident)
+            <|> brackets (sepEndBy (char '$' >> ident) whiteSpace) 
 
-indexForVar :: Parser (Index String)
-indexForVar = try (char '~' >> Superscript <$> ident)
-        <|> try (char '_' >> Subscript <$> ident)
+varNames' :: Parser [Var]
+varNames' = return <$> (char '$' >> identVar)
+            <|> brackets (sepEndBy (char '$' >> identVar) whiteSpace)
 
 argNames :: Parser [Arg]
 argNames = return <$> argName
@@ -500,7 +501,7 @@
     _ | all null vars ->
         let args' = rights args
             args'' = map f (zip args (annonVars 1 (length args)))
-            args''' = map (VarExpr . (either id id)) args''
+            args''' = map (VarExpr . stringToVar . (either id id)) args''
         in return $ ApplyExpr (LambdaExpr (map ScalarArg (rights args'')) (LambdaExpr (map ScalarArg (lefts args'')) $ ApplyExpr func $ TupleExpr args''')) $ TupleExpr args'
       | all (not . null) vars ->
         let ns = Set.fromList $ map read vars
@@ -509,7 +510,7 @@
              then
                let args' = rights args
                    args'' = map g (zip args (annonVars (n + 1) (length args)))
-                   args''' = map (VarExpr . (either id id)) args''
+                   args''' = map (VarExpr . stringToVar . (either id id)) args''
                in return $ ApplyExpr (LambdaExpr (map ScalarArg (rights args'')) (LambdaExpr (map ScalarArg (annonVars 1 n)) $ ApplyExpr func $ TupleExpr args''')) $ TupleExpr args'
              else fail "invalid partial application"
       | otherwise -> fail "invalid partial application"
@@ -581,11 +582,17 @@
 pmapExpr = keywordPmap >> PmapExpr <$> expr <*> expr
 
 subrefsExpr :: Parser EgisonExpr
-subrefsExpr = keywordSubrefs >> SubrefsExpr <$> expr <*> expr
+subrefsExpr = (keywordSubrefs >> SubrefsExpr False <$> expr <*> expr)
+               <|> (keywordSubrefsNew >> SubrefsExpr True <$> expr <*> expr)
 
 suprefsExpr :: Parser EgisonExpr
-suprefsExpr = keywordSuprefs >> SuprefsExpr <$> expr <*> expr
+suprefsExpr = (keywordSuprefs >> SuprefsExpr False <$> expr <*> expr)
+               <|> (keywordSuprefsNew >> SuprefsExpr True <$> expr <*> expr)
 
+userrefsExpr :: Parser EgisonExpr
+userrefsExpr = (keywordUserrefs >> UserrefsExpr False <$> expr <*> expr)
+                <|> (keywordUserrefsNew >> UserrefsExpr True <$> expr <*> expr)
+
 -- Patterns
 
 pattern :: Parser EgisonPattern
@@ -624,7 +631,7 @@
 wildCard = reservedOp "_" >> pure WildCard
 
 patVar :: Parser EgisonPattern
-patVar = PatVar <$> varName
+patVar = char '$' >> PatVar <$> identVarWithoutIndex
 
 varPat :: Parser EgisonPattern
 varPat = VarPat <$> ident
@@ -666,7 +673,7 @@
 dApplyPat = DApplyPat <$> pattern'' <*> sepEndBy pattern whiteSpace 
 
 loopPat :: Parser EgisonPattern
-loopPat = keywordLoop >> LoopPat <$> varName <*> loopRange <*> pattern <*> option (NotPat WildCard) pattern
+loopPat = keywordLoop >> char '$' >> LoopPat <$> identVarWithoutIndex <*> loopRange <*> pattern <*> option (NotPat WildCard) pattern
 
 loopRange :: Parser LoopRange
 loopRange = brackets (try (do s <- expr
@@ -675,7 +682,7 @@
                               return (LoopRange s e ep))
                  <|> (do s <- expr
                          ep <- option WildCard pattern
-                         return (LoopRange s (ApplyExpr (VarExpr "from") (ApplyExpr (VarExpr "-'") (TupleExpr [s, (IntegerExpr 1)]))) ep)))
+                         return (LoopRange s (ApplyExpr (VarExpr $ stringToVar "from") (ApplyExpr (VarExpr $ stringToVar "-'") (TupleExpr [s, (IntegerExpr 1)]))) ep)))
 
 divPat :: Parser EgisonPattern
 divPat = reservedOp "/" >> DivPat <$> pattern <*> pattern
@@ -769,7 +776,7 @@
                 , P.caseSensitive      = True }
 
 symbol0 = oneOf "^"
-symbol1 = oneOf "+-*/=.∂∇"
+symbol1 = oneOf "+-*/.=∂∇"
 symbol2 = symbol1 <|> oneOf "'!?"
 
 lexer :: P.GenTokenParser String () Identity
@@ -823,7 +830,13 @@
   , "pseq"
   , "pmap"
   , "subrefs"
+  , "subrefs!"
   , "suprefs"
+  , "suprefs!"
+  , "user-refs"
+  , "user-refs!"
+  , "function"
+  , "symbolic-tensor"
   , "something"
   , "undefined"]
   
@@ -834,7 +847,6 @@
   , "_"
   , "^"
   , "&"
-  , "|"
   , "|*"
 --  , "'"
 --  , "~"
@@ -904,7 +916,13 @@
 keywordPseq                 = reserved "pseq"
 keywordPmap                 = reserved "pmap"
 keywordSubrefs              = reserved "subrefs"
+keywordSubrefsNew           = reserved "subrefs!"
 keywordSuprefs              = reserved "suprefs"
+keywordSuprefsNew           = reserved "suprefs!"
+keywordUserrefs             = reserved "user-refs"
+keywordUserrefsNew          = reserved "user-refs!"
+keywordFunction             = reserved "function"
+keywordSymbolicTensor       = reserved "symbolic-tensor"
 
 sign :: Num a => Parser (a -> a)
 sign = (char '-' >> return negate)
@@ -961,6 +979,31 @@
 
 ident :: Parser String
 ident = P.identifier lexer
+
+identVar :: Parser Var
+identVar = P.lexeme lexer (do
+  name <- ident
+  is <- many indexType
+  return $ Var (splitOn "." name) is)
+
+identVarWithoutIndex :: Parser Var
+identVarWithoutIndex = do
+    x <- ident
+    return $ stringToVar x
+
+identVarWithIndices :: Parser VarWithIndices
+identVarWithIndices = P.lexeme lexer (do
+  name <- ident
+  is <- many indexForVar
+  return $ VarWithIndices (splitOn "." name) is)
+
+indexForVar :: Parser (Index String)
+indexForVar = try (char '~' >> Superscript <$> ident)
+        <|> try (char '_' >> Subscript <$> ident)
+
+indexType :: Parser (Index ())
+indexType = try (char '~' >> return (Superscript ()))
+        <|> try (char '_' >> return (Subscript ()))
 
 upperName :: Parser String
 upperName = P.lexeme lexer $ upperName'
diff --git a/hs-src/Language/Egison/Primitives.hs b/hs-src/Language/Egison/Primitives.hs
--- a/hs-src/Language/Egison/Primitives.hs
+++ b/hs-src/Language/Egison/Primitives.hs
@@ -45,7 +45,7 @@
   let ops = map (\(name, fn) -> (name, PrimitiveFunc name fn)) (primitives ++ ioPrimitives)
   bindings <- forM (constants ++ ops) $ \(name, op) -> do
     ref <- newIORef . WHNF $ Value op
-    return (name, ref)
+    return (stringToVar name, ref)
   return $ extendEnv nullEnv bindings
 
 primitiveEnvNoIO :: IO Env
@@ -53,64 +53,62 @@
   let ops = map (\(name, fn) -> (name, PrimitiveFunc name fn)) primitives
   bindings <- forM (constants ++ ops) $ \(name, op) -> do
     ref <- newIORef . WHNF $ Value op
-    return (name, ref)
+    return (stringToVar name, ref)
   return $ extendEnv nullEnv bindings
 
 {-# INLINE noArg #-}
 noArg :: EgisonM EgisonValue -> PrimitiveFunc
-noArg f = \args -> do
-  args' <- tupleToList args
-  case args' of 
-    [] -> f >>= return . Value
-    _ -> throwError $ ArgumentsNumPrimitive 0 $ length args'
+noArg f args = do
+    args' <- tupleToList args
+    case args' of
+      [] -> Value <$> f
+      _ -> throwError $ ArgumentsNumPrimitive 0 $ length args'
 
 {-# INLINE oneArg #-}
 oneArg :: (EgisonValue -> EgisonM EgisonValue) -> PrimitiveFunc
-oneArg f = \arg -> do
+oneArg f arg = do
   arg' <- evalWHNF arg
   case arg' of
     (TensorData (Tensor ns ds js)) -> do
-      ds' <- V.mapM (\d -> f d) ds
-      fromTensor (Tensor ns ds' js) >>= return . Value 
-    _ -> f arg' >>= return . Value
+      ds' <- V.mapM f ds
+      Value <$> fromTensor (Tensor ns ds' js)
+    _ -> Value <$> f arg'
 
 {-# INLINE oneArg' #-}
 oneArg' :: (EgisonValue -> EgisonM EgisonValue) -> PrimitiveFunc
-oneArg' f = \arg -> do
+oneArg' f arg = do
   arg' <- evalWHNF arg
-  case arg' of
-    _ -> f arg' >>= return . Value
+  Value <$> f arg'
 
 {-# INLINE twoArgs #-}
 twoArgs :: (EgisonValue -> EgisonValue -> EgisonM EgisonValue) -> PrimitiveFunc
-twoArgs f = \args -> do
+twoArgs f args = do
   args' <- tupleToList args
   case args' of 
-    [(TensorData t1@(Tensor _ _ _)), (TensorData t2@(Tensor _ _ _))] -> do
-      tProduct f t1 t2 >>= fromTensor >>= return . Value
-    [(TensorData(Tensor ns ds js)), val] -> do
+    [TensorData t1@(Tensor _ _ _), TensorData t2@(Tensor _ _ _)] -> Value <$> (tProduct f t1 t2 >>= fromTensor)
+    [TensorData(Tensor ns ds js), val] -> do
       ds' <- V.mapM (\d -> f d val) ds
-      fromTensor (Tensor ns ds' js) >>= return . Value 
-    [val, (TensorData (Tensor ns ds js))] -> do
+      Value <$> fromTensor (Tensor ns ds' js)
+    [val, TensorData (Tensor ns ds js)] -> do
       ds' <- V.mapM (\d -> f val d) ds
-      fromTensor (Tensor ns ds' js) >>= return . Value 
-    [val, val'] -> f val val' >>= return . Value
+      Value <$> fromTensor (Tensor ns ds' js)
+    [val, val'] -> Value <$> f val val'
     _ -> throwError $ ArgumentsNumPrimitive 2 $ length args'
 
 {-# INLINE twoArgs' #-}
 twoArgs' :: (EgisonValue -> EgisonValue -> EgisonM EgisonValue) -> PrimitiveFunc
-twoArgs' f = \args -> do
+twoArgs' f args = do
   args' <- tupleToList args
   case args' of 
-    [val, val'] -> f val val' >>= return . Value
+    [val, val'] -> Value <$> f val val'
     _ -> throwError $ ArgumentsNumPrimitive 2 $ length args'
 
 {-# INLINE threeArgs' #-}
 threeArgs' :: (EgisonValue -> EgisonValue -> EgisonValue -> EgisonM EgisonValue) -> PrimitiveFunc
-threeArgs' f = \args -> do
+threeArgs' f args = do
   args' <- tupleToList args
   case args' of 
-    [val, val', val''] -> f val val' val'' >>= return . Value
+    [val, val', val''] -> Value <$> f val val' val''
     _ -> throwError $ ArgumentsNumPrimitive 3 $ length args'
 
 --
@@ -145,8 +143,6 @@
              , ("from-math-expr", fromScalarData)
              , ("to-math-expr", toScalarData)
              , ("to-math-expr'", toScalarData)
-             , ("decons-user-scripts", deconsUserScripts)
-             , ("append-user-scripts", appendUserScripts)
 
              , ("modulo",    integerBinaryOp mod)
              , ("quotient",   integerBinaryOp quot)
@@ -244,7 +240,7 @@
 rationalBinaryOp op = twoArgs $ \val val' -> do
   r <- fromEgison val :: EgisonM Rational
   r' <- fromEgison val' :: EgisonM Rational
-  let r'' = (op r r'')
+  let r'' = op r r''
   return $ toEgison r''
 
 rationalBinaryPred :: (Rational -> Rational -> Bool) -> PrimitiveFunc
@@ -266,16 +262,14 @@
   return $ Bool $ pred i i'
 
 floatUnaryOp :: (Double -> Double) -> PrimitiveFunc
-floatUnaryOp op = oneArg $ \val -> do
-  case val of
-    (Float f 0) -> return $ Float (op f) 0
-    _ -> throwError $ TypeMismatch "float" (Value val)
+floatUnaryOp op = oneArg $ \val -> case val of
+                                     (Float f 0) -> return $ Float (op f) 0
+                                     _ -> throwError $ TypeMismatch "float" (Value val)
 
 floatBinaryOp :: (Double -> Double -> Double) -> PrimitiveFunc
-floatBinaryOp op = twoArgs $ \val val' -> do
-  case (val, val') of
-    ((Float f 0), (Float f' 0)) -> return $ Float (op f f') 0
-    _ -> throwError $ TypeMismatch "float" (Value val)
+floatBinaryOp op = twoArgs $ \val val' -> case (val, val') of
+                                            (Float f 0, Float f' 0) -> return $ Float (op f f') 0
+                                            _ -> throwError $ TypeMismatch "float" (Value val)
 
 floatBinaryPred :: (Double -> Double -> Bool) -> PrimitiveFunc
 floatBinaryPred pred = twoArgs $ \val val' -> do
@@ -284,28 +278,24 @@
   return $ Bool $ pred f f'
 
 floatPlus :: PrimitiveFunc
-floatPlus = twoArgs $ \val val' -> do
-  case (val, val') of
-    ((Float x y), (Float x' y')) -> return $ Float (x + x')  (y + y')
-    _ -> throwError $ TypeMismatch "float" (Value val)
+floatPlus = twoArgs $ \val val' -> case (val, val') of
+                                     (Float x y, Float x' y') -> return $ Float (x + x')  (y + y')
+                                     _ -> throwError $ TypeMismatch "float" (Value val)
 
 floatMinus :: PrimitiveFunc
-floatMinus = twoArgs $ \val val' -> do
-  case (val, val') of
-    ((Float x y), (Float x' y')) -> return $ Float (x - x')  (y - y')
-    _ -> throwError $ TypeMismatch "float" (Value val)
+floatMinus = twoArgs $ \val val' -> case (val, val') of
+                                      (Float x y, Float x' y') -> return $ Float (x - x')  (y - y')
+                                      _ -> throwError $ TypeMismatch "float" (Value val)
 
 floatMult :: PrimitiveFunc
-floatMult = twoArgs $ \val val' -> do
-  case (val, val') of
-    ((Float x y), (Float x' y')) -> return $ Float (x * x' - y * y')  (x * y' + x' * y)
-    _ -> throwError $ TypeMismatch "float" (Value val)
+floatMult = twoArgs $ \val val' -> case (val, val') of
+                                     (Float x y, Float x' y') -> return $ Float (x * x' - y * y')  (x * y' + x' * y)
+                                     _ -> throwError $ TypeMismatch "float" (Value val)
 
 floatDivide :: PrimitiveFunc
-floatDivide = twoArgs $ \val val' -> do
-  case (val, val') of
-    ((Float x y), (Float x' y')) -> return $ Float ((x * x' + y * y') / (x' * x' + y' * y')) ((y * x' - x * y') / (x' * x' + y' * y'))
-    _ -> throwError $ TypeMismatch "float" (Value val)
+floatDivide = twoArgs $ \val val' -> case (val, val') of
+                                       (Float x y, Float x' y') -> return $ Float ((x * x' + y * y') / (x' * x' + y' * y')) ((y * x' - x * y') / (x' * x' + y' * y'))
+                                       _ -> throwError $ TypeMismatch "float" (Value val)
 
 
 --
@@ -313,10 +303,10 @@
 --
 
 scalarBinaryOp :: (ScalarData -> ScalarData -> ScalarData) -> PrimitiveFunc
-scalarBinaryOp mOp = twoArgs $ \val val' -> do
-  scalarBinaryOp' val val'
+scalarBinaryOp mOp = twoArgs $ \val val' -> scalarBinaryOp' val val'
  where
   scalarBinaryOp' (ScalarData m1) (ScalarData m2) = (return . ScalarData . mathNormalize') (mOp m1 m2)
+  scalarBinaryOp' (ScalarData _)  val             = throwError $ TypeMismatch "number" (Value val)
   scalarBinaryOp' val             _               = throwError $ TypeMismatch "number" (Value val)
 
 plus :: PrimitiveFunc
@@ -332,49 +322,37 @@
 divide = scalarBinaryOp (\m1 (Div p1 p2) -> mathMult m1 (Div p2 p1))
 
 numerator' :: PrimitiveFunc
-numerator' =  oneArg $ numerator''
+numerator' =  oneArg numerator''
  where
   numerator'' (ScalarData m) = return $ ScalarData (mathNumerator m)
   numerator'' val = throwError $ TypeMismatch "rational" (Value val)
 
 denominator' :: PrimitiveFunc
-denominator' =  oneArg $ denominator''
+denominator' =  oneArg denominator''
  where
   denominator'' (ScalarData m) = return $ ScalarData (mathDenominator m)
   denominator'' val = throwError $ TypeMismatch "rational" (Value val)
 
 fromScalarData :: PrimitiveFunc
-fromScalarData = oneArg $ fromScalarData'
+fromScalarData = oneArg fromScalarData'
  where
   fromScalarData' (ScalarData m) = return $ mathExprToEgison m
   fromScalarData' val = throwError $ TypeMismatch "number" (Value val)
 
 toScalarData :: PrimitiveFunc
-toScalarData = oneArg $ toScalarData'
- where
-  toScalarData' val = egisonToScalarData val >>= return . ScalarData . mathNormalize'
-
-appendUserScripts :: PrimitiveFunc
-appendUserScripts = twoArgs $ appendUserScripts'
+toScalarData = oneArg toScalarData'
  where
-  appendUserScripts' v (Collection is) = do
-    let is' = map Userscript (toList is)
-    return $ UserIndexedData v is'
+  toScalarData' val = (ScalarData . mathNormalize') <$> egisonToScalarData val
 
-deconsUserScripts :: PrimitiveFunc
-deconsUserScripts = oneArg $ deconsUserScripts'
- where
-  deconsUserScripts' (UserIndexedData v is) = return $ Tuple [v, Collection (Sq.fromList (map (\(Userscript i) -> i) is))]
-  deconsUserScripts' v = return $ Tuple [v, Collection (Sq.fromList [])]
 --
 -- Pred
 --
 eq :: PrimitiveFunc
-eq = twoArgs $ \val val' ->
+eq = twoArgs' $ \val val' ->
   return $ Bool $ val == val'
 
 lt :: PrimitiveFunc
-lt = twoArgs $ \val val' -> scalarBinaryPred' val val'
+lt = twoArgs' $ \val val' -> scalarBinaryPred' val val'
  where
   scalarBinaryPred' m@(ScalarData _) n@(ScalarData _) = do
     r <- fromEgison m :: EgisonM Rational
@@ -386,7 +364,7 @@
   scalarBinaryPred' val          _             = throwError $ TypeMismatch "number" (Value val)
   
 lte :: PrimitiveFunc
-lte = twoArgs $ \val val' -> scalarBinaryPred' val val'
+lte = twoArgs' $ \val val' -> scalarBinaryPred' val val'
  where
   scalarBinaryPred' m@(ScalarData _) n@(ScalarData _) = do
     r <- fromEgison m :: EgisonM Rational
@@ -398,7 +376,7 @@
   scalarBinaryPred' val          _             = throwError $ TypeMismatch "number" (Value val)
   
 gt :: PrimitiveFunc
-gt = twoArgs $ \val val' -> scalarBinaryPred' val val'
+gt = twoArgs' $ \val val' -> scalarBinaryPred' val val'
  where
   scalarBinaryPred' m@(ScalarData _) n@(ScalarData _) = do
     r <- fromEgison m :: EgisonM Rational
@@ -410,7 +388,7 @@
   scalarBinaryPred' val          _             = throwError $ TypeMismatch "number" (Value val)
   
 gte :: PrimitiveFunc
-gte = twoArgs $ \val val' -> scalarBinaryPred' val val'
+gte = twoArgs' $ \val val' -> scalarBinaryPred' val val'
  where
   scalarBinaryPred' m@(ScalarData _) n@(ScalarData _) = do
     r <- fromEgison m :: EgisonM Rational
@@ -425,18 +403,18 @@
 truncate' = oneArg $ \val -> numberUnaryOp' val
  where
   numberUnaryOp' (ScalarData (Div (Plus []) _)) = return $ toEgison (0 :: Integer)
-  numberUnaryOp' (ScalarData (Div (Plus [(Term x [])]) (Plus [(Term y [])]))) = return $ toEgison (quot x y)
-  numberUnaryOp' (Float x _)           = return $ toEgison ((truncate x) :: Integer)
+  numberUnaryOp' (ScalarData (Div (Plus [Term x []]) (Plus [Term y []]))) = return $ toEgison (quot x y)
+  numberUnaryOp' (Float x _)           = return $ toEgison (truncate x :: Integer)
   numberUnaryOp' val                   = throwError $ TypeMismatch "rational or float" (Value val)
 
 realPart :: PrimitiveFunc
-realPart =  oneArg $ realPart'
+realPart =  oneArg realPart'
  where
   realPart' (Float x y) = return $ Float x 0
   realPart' val = throwError $ TypeMismatch "float" (Value val)
 
 imaginaryPart :: PrimitiveFunc
-imaginaryPart =  oneArg $ imaginaryPart'
+imaginaryPart =  oneArg imaginaryPart'
  where
   imaginaryPart' (Float _ y) = return $ Float y 0
   imaginaryPart' val = throwError $ TypeMismatch "float" (Value val)
@@ -446,21 +424,21 @@
 --
 
 tensorSize' :: PrimitiveFunc
-tensorSize' = oneArg' $ tensorSize''
+tensorSize' = oneArg' tensorSize''
  where
   tensorSize'' (TensorData (Tensor ns _ _)) = return . Collection . Sq.fromList $ map toEgison ns
-  tensorSize'' _ = return . Collection $ Sq.fromList $ []
+  tensorSize'' _ = return . Collection $ Sq.fromList []
 
 tensorToList' :: PrimitiveFunc
-tensorToList' = oneArg' $ tensorToList''
+tensorToList' = oneArg' tensorToList''
  where
   tensorToList'' (TensorData (Tensor _ xs _)) = return . Collection . Sq.fromList $ V.toList xs
-  tensorToList'' x = return . Collection $ Sq.fromList $ [x]
+  tensorToList'' x = return . Collection $ Sq.fromList [x]
 
 dfOrder' :: PrimitiveFunc
-dfOrder' = oneArg' $ dfOrder''
+dfOrder' = oneArg' dfOrder''
  where
-  dfOrder'' (TensorData (Tensor ns _ is)) = return (toEgison ((fromIntegral ((length ns) - (length is))) :: Integer))
+  dfOrder'' (TensorData (Tensor ns _ is)) = return (toEgison (fromIntegral (length ns - length is) :: Integer))
   dfOrder'' _ = return (toEgison (0 :: Integer))
 
 --
@@ -468,7 +446,7 @@
 --
 numberToFloat' :: EgisonValue -> EgisonValue
 numberToFloat' (ScalarData (Div (Plus []) _)) = Float 0 0
-numberToFloat' (ScalarData (Div (Plus [(Term x [])]) (Plus [(Term y [])]))) = Float (fromRational (x % y)) 0
+numberToFloat' (ScalarData (Div (Plus [Term x []]) (Plus [Term y []]))) = Float (fromRational (x % y)) 0
 
 integerToFloat :: PrimitiveFunc
 integerToFloat = rationalToFloat
@@ -477,24 +455,22 @@
 rationalToFloat = oneArg $ \val ->
   case val of
     (ScalarData (Div (Plus []) _)) -> return $ numberToFloat' val
-    (ScalarData (Div (Plus [(Term _ [])]) (Plus [(Term _ [])]))) -> return $ numberToFloat' val
+    (ScalarData (Div (Plus [Term _ []]) (Plus [Term _ []]))) -> return $ numberToFloat' val
     _ -> throwError $ TypeMismatch "integer or rational number" (Value val)
 
 charToInteger :: PrimitiveFunc
-charToInteger = oneArg $ \val -> do
-  case val of
-    Char c -> do
-      let i = fromIntegral $ ord c :: Integer
-      return $ toEgison i
-    _ -> throwError $ TypeMismatch "character" (Value val)
+charToInteger = oneArg $ \val -> case val of
+                                   Char c -> do
+                                     let i = fromIntegral $ ord c :: Integer
+                                     return $ toEgison i
+                                   _ -> throwError $ TypeMismatch "character" (Value val)
 
 integerToChar :: PrimitiveFunc
-integerToChar = oneArg $ \val -> do
-  case val of
-    (ScalarData _) -> do
-       i <- fromEgison val :: EgisonM Integer
-       return $ Char $ chr $ fromIntegral i
-    _ -> throwError $ TypeMismatch "integer" (Value val)
+integerToChar = oneArg $ \val -> case val of
+                                   (ScalarData _) -> do
+                                      i <- fromEgison val :: EgisonM Integer
+                                      return $ Char $ chr $ fromIntegral i
+                                   _ -> throwError $ TypeMismatch "integer" (Value val)
 
 floatToIntegerOp :: (Double -> Integer) -> PrimitiveFunc
 floatToIntegerOp op = oneArg $ \val -> do
@@ -510,61 +486,54 @@
   return $ String str
 
 unpack :: PrimitiveFunc
-unpack = oneArg $ \val -> do
-  case val of
-    String str -> return $ toEgison (T.unpack str)
-    _ -> throwError $ TypeMismatch "string" (Value val)
+unpack = oneArg $ \val -> case val of
+                            String str -> return $ toEgison (T.unpack str)
+                            _ -> throwError $ TypeMismatch "string" (Value val)
 
 unconsString :: PrimitiveFunc
-unconsString = oneArg $ \val -> do
-  case val of
-    String str -> case T.uncons str of
-                    Just (c, rest) ->  return $ Tuple [Char c, String rest]
-                    Nothing -> throwError $ Default "Tried to unsnoc empty string"
-    _ -> throwError $ TypeMismatch "string" (Value val)
+unconsString = oneArg $ \val -> case val of
+                                  String str -> case T.uncons str of
+                                                  Just (c, rest) ->  return $ Tuple [Char c, String rest]
+                                                  Nothing -> throwError $ Default "Tried to unsnoc empty string"
+                                  _ -> throwError $ TypeMismatch "string" (Value val)
 
 lengthString :: PrimitiveFunc
-lengthString = oneArg $ \val -> do
-  case val of
-    String str -> return . (\x -> toEgison x) . toInteger $ T.length str
-    _ -> throwError $ TypeMismatch "string" (Value val)
+lengthString = oneArg $ \val -> case val of
+                                  String str -> return . toEgison . toInteger $ T.length str
+                                  _ -> throwError $ TypeMismatch "string" (Value val)
 
 appendString :: PrimitiveFunc
-appendString = twoArgs $ \val1 val2 -> do
-  case (val1, val2) of
-    (String str1, String str2) -> return . String $ T.append str1 str2
-    (String _, _) -> throwError $ TypeMismatch "string" (Value val2)
-    (_, _) -> throwError $ TypeMismatch "string" (Value val1)
+appendString = twoArgs $ \val1 val2 -> case (val1, val2) of
+                                         (String str1, String str2) -> return . String $ T.append str1 str2
+                                         (String _, _) -> throwError $ TypeMismatch "string" (Value val2)
+                                         (_, _) -> throwError $ TypeMismatch "string" (Value val1)
 
 splitString :: PrimitiveFunc
-splitString = twoArgs $ \pat src -> do
-  case (pat, src) of
-    (String patStr, String srcStr) -> return . Collection . Sq.fromList $ map String $ T.splitOn patStr srcStr
-    (String _, _) -> throwError $ TypeMismatch "string" (Value src)
-    (_, _) -> throwError $ TypeMismatch "string" (Value pat)
+splitString = twoArgs $ \pat src -> case (pat, src) of
+                                      (String patStr, String srcStr) -> return . Collection . Sq.fromList $ map String $ T.splitOn patStr srcStr
+                                      (String _, _) -> throwError $ TypeMismatch "string" (Value src)
+                                      (_, _) -> throwError $ TypeMismatch "string" (Value pat)
 
 regexString :: PrimitiveFunc
-regexString = twoArgs $ \pat src -> do
-  case (pat, src) of
-    (String patStr, String srcStr) -> do
-      let (a, b, c) = (((T.unpack srcStr) =~ (T.unpack patStr)) :: (String, String, String))
-      if b == ""
-        then return . Collection . Sq.fromList $ []
-        else return . Collection . Sq.fromList $ [Tuple [String $ T.pack a, String $ T.pack b, String $ T.pack c]]
-    (String _, _) -> throwError $ TypeMismatch "string" (Value src)
-    (_, _) -> throwError $ TypeMismatch "string" (Value pat)
+regexString = twoArgs $ \pat src -> case (pat, src) of
+                                      (String patStr, String srcStr) -> do
+                                        let (a, b, c) = (T.unpack srcStr =~ T.unpack patStr) :: (String, String, String)
+                                        if b == ""
+                                          then return . Collection . Sq.fromList $ []
+                                          else return . Collection . Sq.fromList $ [Tuple [String $ T.pack a, String $ T.pack b, String $ T.pack c]]
+                                      (String _, _) -> throwError $ TypeMismatch "string" (Value src)
+                                      (_, _) -> throwError $ TypeMismatch "string" (Value pat)
 
 regexStringCaptureGroup :: PrimitiveFunc
-regexStringCaptureGroup = twoArgs $ \pat src -> do
-  case (pat, src) of
-    (String patStr, String srcStr) -> do
-      let ret = (((T.unpack srcStr) =~ (T.unpack patStr)) :: [[String]])
-      case ret of 
-        [] -> return . Collection . Sq.fromList $ []
-        ((x:xs):_) -> do let (a, c) = T.breakOn (T.pack x) srcStr
-                         return . Collection . Sq.fromList $ [Tuple [String a, Collection (Sq.fromList (map (String . T.pack) xs)), String (T.drop (length x) c)]]
-    (String _, _) -> throwError $ TypeMismatch "string" (Value src)
-    (_, _) -> throwError $ TypeMismatch "string" (Value pat)
+regexStringCaptureGroup = twoArgs $ \pat src -> case (pat, src) of
+                                                  (String patStr, String srcStr) -> do
+                                                    let ret = (T.unpack srcStr =~ T.unpack patStr) :: [[String]]
+                                                    case ret of 
+                                                      [] -> return . Collection . Sq.fromList $ []
+                                                      ((x:xs):_) -> do let (a, c) = T.breakOn (T.pack x) srcStr
+                                                                       return . Collection . Sq.fromList $ [Tuple [String a, Collection (Sq.fromList (map (String . T.pack) xs)), String (T.drop (length x) c)]]
+                                                  (String _, _) -> throwError $ TypeMismatch "string" (Value src)
+                                                  (_, _) -> throwError $ TypeMismatch "string" (Value pat)
 
 --regexStringMatch :: PrimitiveFunc
 --regexStringMatch = twoArgs $ \pat src -> do
@@ -574,42 +543,38 @@
 --    (_, _) -> throwError $ TypeMismatch "string" (Value pat)
 
 addPrime :: PrimitiveFunc
-addPrime = oneArg $ \sym -> do
-  case sym of
-    ScalarData (Div (Plus [(Term 1 [(Symbol id name is, 1)])]) (Plus [(Term 1 [])])) -> return (ScalarData (Div (Plus [(Term 1 [(Symbol id (name ++ "'") is, 1)])]) (Plus [(Term 1 [])])))
-    _ ->  throwError $ TypeMismatch "symbol" (Value sym)
+addPrime = oneArg $ \sym -> case sym of
+                              ScalarData (Div (Plus [Term 1 [(Symbol id name is, 1)]]) (Plus [Term 1 []])) -> return (ScalarData (Div (Plus [Term 1 [(Symbol id (name ++ "'") is, 1)]]) (Plus [Term 1 []])))
+                              _ ->  throwError $ TypeMismatch "symbol" (Value sym)
 
 addSubscript :: PrimitiveFunc
-addSubscript = twoArgs $ \fn sub -> do
-  case (fn, sub) of
-    (ScalarData (Div (Plus [(Term 1 [(Symbol id name is, 1)])]) (Plus [(Term 1 [])])),
-     ScalarData s@(Div (Plus [(Term 1 [(Symbol _ _ [], 1)])]) (Plus [(Term 1 [])]))) -> return (ScalarData (Div (Plus [(Term 1 [(Symbol id name (is ++ [Subscript s]), 1)])]) (Plus [(Term 1 [])])))
-    (ScalarData (Div (Plus [(Term 1 [(Symbol id name is, 1)])]) (Plus [(Term 1 [])])),
-     ScalarData s@(Div (Plus [(Term _ [])]) (Plus [(Term 1 [])]))) -> return (ScalarData (Div (Plus [(Term 1 [(Symbol id name (is ++ [Subscript s]), 1)])]) (Plus [(Term 1 [])])))
-    (ScalarData (Div (Plus [(Term 1 [(Symbol _ _ _, 1)])]) (Plus [(Term 1 [])])),
-     _) -> throwError $ TypeMismatch "symbol or integer" (Value sub)
-    _ ->  throwError $ TypeMismatch "symbol or integer" (Value fn)
+addSubscript = twoArgs $ \fn sub -> case (fn, sub) of
+                                      (ScalarData (Div (Plus [Term 1 [(Symbol id name is, 1)]]) (Plus [Term 1 []])),
+                                       ScalarData s@(Div (Plus [Term 1 [(Symbol _ _ [], 1)]]) (Plus [Term 1 []]))) -> return (ScalarData (Div (Plus [Term 1 [(Symbol id name (is ++ [Subscript s]), 1)]]) (Plus [Term 1 []])))
+                                      (ScalarData (Div (Plus [Term 1 [(Symbol id name is, 1)]]) (Plus [Term 1 []])),
+                                       ScalarData s@(Div (Plus [Term _ []]) (Plus [Term 1 []]))) -> return (ScalarData (Div (Plus [Term 1 [(Symbol id name (is ++ [Subscript s]), 1)]]) (Plus [Term 1 []])))
+                                      (ScalarData (Div (Plus [Term 1 [(Symbol _ _ _, 1)]]) (Plus [Term 1 []])),
+                                       _) -> throwError $ TypeMismatch "symbol or integer" (Value sub)
+                                      _ ->  throwError $ TypeMismatch "symbol or integer" (Value fn)
 
 addSuperscript :: PrimitiveFunc
-addSuperscript = twoArgs $ \fn sub -> do
-  case (fn, sub) of
-    (ScalarData (Div (Plus [(Term 1 [(Symbol id name is, 1)])]) (Plus [(Term 1 [])])),
-     ScalarData s@(Div (Plus [(Term 1 [(Symbol _ _ [], 1)])]) (Plus [(Term 1 [])]))) -> return (ScalarData (Div (Plus [(Term 1 [(Symbol id name (is ++ [Superscript s]), 1)])]) (Plus [(Term 1 [])])))
-    (ScalarData (Div (Plus [(Term 1 [(Symbol id name is, 1)])]) (Plus [(Term 1 [])])),
-     ScalarData s@(Div (Plus [(Term _ [])]) (Plus [(Term 1 [])]))) -> return (ScalarData (Div (Plus [(Term 1 [(Symbol id name (is ++ [Superscript s]), 1)])]) (Plus [(Term 1 [])])))
-    (ScalarData (Div (Plus [(Term 1 [(Symbol _ _ _, 1)])]) (Plus [(Term 1 [])])),
-     _) -> throwError $ TypeMismatch "symbol" (Value sub)
-    _ ->  throwError $ TypeMismatch "symbol" (Value fn)
+addSuperscript = twoArgs $ \fn sub -> case (fn, sub) of
+                                        (ScalarData (Div (Plus [Term 1 [(Symbol id name is, 1)]]) (Plus [Term 1 []])),
+                                         ScalarData s@(Div (Plus [Term 1 [(Symbol _ _ [], 1)]]) (Plus [Term 1 []]))) -> return (ScalarData (Div (Plus [Term 1 [(Symbol id name (is ++ [Superscript s]), 1)]]) (Plus [Term 1 []])))
+                                        (ScalarData (Div (Plus [Term 1 [(Symbol id name is, 1)]]) (Plus [Term 1 []])),
+                                         ScalarData s@(Div (Plus [Term _ []]) (Plus [Term 1 []]))) -> return (ScalarData (Div (Plus [Term 1 [(Symbol id name (is ++ [Superscript s]), 1)]]) (Plus [Term 1 []])))
+                                        (ScalarData (Div (Plus [Term 1 [(Symbol _ _ _, 1)]]) (Plus [Term 1 []])),
+                                         _) -> throwError $ TypeMismatch "symbol" (Value sub)
+                                        _ ->  throwError $ TypeMismatch "symbol" (Value fn)
 
 readProcess' :: PrimitiveFunc
-readProcess' = threeArgs' $ \cmd args input -> do
-  case (cmd, args, input) of
-    (String cmdStr, Collection argStrs, String inputStr) -> do
-      outputStr <- liftIO $ readProcess (T.unpack cmdStr) (map (\arg -> case arg of
-                                                                          String argStr -> T.unpack argStr)
-                                                                (toList argStrs)) (T.unpack inputStr)
-      return (String (T.pack outputStr))
-    (_, _, _) -> throwError $ TypeMismatch "(string, collection, string)" (Value (Tuple [cmd, args, input]))
+readProcess' = threeArgs' $ \cmd args input -> case (cmd, args, input) of
+                                                 (String cmdStr, Collection argStrs, String inputStr) -> do
+                                                   outputStr <- liftIO $ readProcess (T.unpack cmdStr) (map (\arg -> case arg of
+                                                                                                                       String argStr -> T.unpack argStr)
+                                                                                                             (toList argStrs)) (T.unpack inputStr)
+                                                   return (String (T.pack outputStr))
+                                                 (_, _, _) -> throwError $ TypeMismatch "(string, collection, string)" (Value (Tuple [cmd, args, input]))
 
 read' :: PrimitiveFunc
 read'= oneArg' $ \val -> fromEgison val >>= readExpr . T.unpack >>= evalExprDeep nullEnv
@@ -641,14 +606,14 @@
   mRet <- runMaybeT (unconsCollection whnf)
   case mRet of
     Just (carObjRef, cdrObjRef) -> return $ Intermediate $ ITuple [carObjRef, cdrObjRef]
-    Nothing -> throwError $ Default $ "cannot uncons collection"
+    Nothing -> throwError $ Default "cannot uncons collection"
 
 unsnoc' :: PrimitiveFunc
 unsnoc' whnf = do
   mRet <- runMaybeT (unsnocCollection whnf)
   case mRet of
     Just (racObjRef, rdcObjRef) -> return $ Intermediate $ ITuple [racObjRef, rdcObjRef]
-    Nothing -> throwError $ Default $ "cannot unsnoc collection"
+    Nothing -> throwError $ Default "cannot unsnoc collection"
 
 -- Test
 
@@ -660,11 +625,10 @@
     else throwError $ Assertion $ show label
 
 assertEqual :: PrimitiveFunc
-assertEqual = threeArgs' $ \label actual expected -> do
-  if actual == expected
-    then return $ Bool True
-    else throwError $ Assertion $ show label ++ "\n expected: " ++ show expected ++
-                                  "\n but found: " ++ show actual
+assertEqual = threeArgs' $ \label actual expected -> if actual == expected
+                                                       then return $ Bool True
+                                                       else throwError $ Assertion $ show label ++ "\n expected: " ++ show expected ++
+                                                                                      "\n but found: " ++ show actual
 
 --
 -- IO Primitives
diff --git a/hs-src/Language/Egison/Types.hs b/hs-src/Language/Egison/Types.hs
--- a/hs-src/Language/Egison/Types.hs
+++ b/hs-src/Language/Egison/Types.hs
@@ -1,6 +1,6 @@
 {-# Language TypeSynonymInstances, FlexibleInstances, GeneralizedNewtypeDeriving,
              MultiParamTypeClasses, UndecidableInstances, DeriveDataTypeable,
-             TypeFamilies, TupleSections #-}
+             TypeFamilies, TupleSections, DeriveGeneric #-}
 {- |
 Module      : Language.Egison.Types
 Copyright   : Satoshi Egi
@@ -17,7 +17,6 @@
     , EgisonPattern (..)
     , Arg (..)
     , Index (..)
-    , UserIndex (..)
     , InnerExpr (..)
     , BindingExpr (..)
     , MatchClause (..)
@@ -60,6 +59,8 @@
     , tConcat'
     -- * Scalar
     , symbolScalarData
+    , getSymId
+    , getSymName
     , mathExprToEgison
     , egisonToScalarData
     , mathNormalize'
@@ -145,6 +146,8 @@
     , isArray'
     , isHash'
     , readUTF8File
+    , stringToVar
+    , varToVarWithIndices
     ) where
 
 import Prelude hiding (foldr, mappend, mconcat)
@@ -169,11 +172,13 @@
 import Data.Sequence (Seq)
 import Data.Foldable (foldr, toList)
 import Data.IORef
+import Data.Hashable (Hashable)
 import Data.HashMap.Strict (HashMap)
 import qualified Data.HashMap.Strict as HashMap
 
-import Data.List (intercalate, sort, sortBy, find, findIndex, splitAt, (\\), elem, delete, deleteBy, any, partition)
-import Data.Text (Text)
+import Data.List (intercalate, sort, sortBy, find, findIndex, splitAt, (\\), elem, delete, deleteBy, any, partition, intercalate, elemIndex)
+import Data.List.Split (splitOn)
+import Data.Text (Text, pack)
 import qualified Data.Text as T
 
 import System.IO
@@ -182,6 +187,8 @@
 
 import System.IO.Unsafe (unsafePerformIO)
 
+import GHC.Generics (Generic)
+
 --
 -- Expressions
 --
@@ -202,12 +209,12 @@
   | BoolExpr Bool
   | IntegerExpr Integer
   | FloatExpr Double Double
-  | VarExpr String
+  | VarExpr Var
   | FreshVarExpr
-  | IndexedExpr Bool EgisonExpr [Index EgisonExpr]
-  | SubrefsExpr EgisonExpr EgisonExpr
-  | SuprefsExpr EgisonExpr EgisonExpr
-  | UserIndexedExpr EgisonExpr [UserIndex EgisonExpr]
+  | IndexedExpr Bool EgisonExpr [Index EgisonExpr]  -- True -> delete old index and append new one
+  | SubrefsExpr Bool EgisonExpr EgisonExpr
+  | SuprefsExpr Bool EgisonExpr EgisonExpr
+  | UserrefsExpr Bool EgisonExpr EgisonExpr
   | PowerExpr EgisonExpr EgisonExpr
   | InductiveDataExpr String [EgisonExpr]
   | TupleExpr [EgisonExpr]
@@ -245,7 +252,7 @@
   | AlgebraicDataMatcherExpr [(String, [EgisonExpr])]
 
   | QuoteExpr EgisonExpr
-  | QuoteFunctionExpr EgisonExpr
+  | QuoteSymbolExpr EgisonExpr
   
   | WedgeExpr EgisonExpr
   | WedgeApplyExpr EgisonExpr EgisonExpr
@@ -276,6 +283,9 @@
   | TransposeExpr EgisonExpr EgisonExpr
   | FlipIndicesExpr EgisonExpr
 
+  | FunctionExpr [EgisonExpr]
+  | SymbolicTensorExpr [EgisonExpr] EgisonExpr String
+
   | SomethingExpr
   | UndefinedExpr
  deriving (Eq)
@@ -290,24 +300,24 @@
     Subscript a
   | Superscript a
   | SupSubscript a
+  | MultiSubscript a a
+  | MultiSuperscript a a
   | DFscript Integer Integer -- DifferentialForm
- deriving (Eq)
-
-data UserIndex a = Userscript a
- deriving (Eq)
+  | Userscript a
+ deriving (Eq, Generic)
 
 data InnerExpr =
     ElementExpr EgisonExpr
   | SubCollectionExpr EgisonExpr
  deriving (Show, Eq)
 
-type BindingExpr = ([String], EgisonExpr)
+type BindingExpr = ([Var], EgisonExpr)
 type MatchClause = (EgisonPattern, EgisonExpr)
 type MatcherInfo = [(PrimitivePatPattern, EgisonExpr, [(PrimitiveDataPattern, EgisonExpr)])]
 
 data EgisonPattern =
     WildCard
-  | PatVar String
+  | PatVar Var
   | ValuePat EgisonExpr
   | PredPat EgisonExpr
   | IndexedPat EgisonPattern [EgisonExpr]
@@ -318,7 +328,7 @@
   | OrderedOrPat [EgisonPattern]
   | TuplePat [EgisonPattern]
   | InductivePat String [EgisonPattern]
-  | LoopPat String LoopRange EgisonPattern EgisonPattern
+  | LoopPat Var LoopRange EgisonPattern EgisonPattern
   | ContPat
   | PApplyPat EgisonExpr [EgisonPattern]
   | VarPat String
@@ -362,7 +372,6 @@
   | Bool Bool
   | ScalarData ScalarData
   | TensorData (Tensor EgisonValue)
-  | UserIndexedData EgisonValue [UserIndex EgisonValue]
   | Float Double Double
   | InductiveData String [EgisonValue]
   | Tuple [EgisonValue]
@@ -372,10 +381,10 @@
   | CharHash (HashMap Char EgisonValue)
   | StrHash (HashMap Text EgisonValue)
   | UserMatcher Env PMMode MatcherInfo
-  | Func (Maybe String) Env [String] EgisonExpr
+  | Func (Maybe Var) Env [String] EgisonExpr
   | PartialFunc Env Integer EgisonExpr
-  | CFunc (Maybe String) Env String EgisonExpr
-  | MemoizedFunc (Maybe String) ObjectRef (IORef (HashMap [Integer] ObjectRef)) Env [String] EgisonExpr
+  | CFunc (Maybe Var) Env String EgisonExpr
+  | MemoizedFunc (Maybe Var) ObjectRef (IORef (HashMap [Integer] ObjectRef)) Env [String] EgisonExpr
   | Proc (Maybe String) Env [String] EgisonExpr
   | Macro [String] EgisonExpr
   | PatternFunc Env [String] EgisonPattern
@@ -395,7 +404,7 @@
     Div PolyExpr PolyExpr
  deriving (Eq)
 
-data PolyExpr =
+newtype PolyExpr =
     Plus [TermExpr]
 
 data TermExpr =
@@ -405,14 +414,15 @@
     Symbol String String [Index ScalarData] -- ID, Name, Indices
   | Apply EgisonValue [ScalarData]
   | Quote ScalarData
+  | FunctionData (Maybe EgisonValue) [EgisonValue] [EgisonValue] [Index ScalarData] -- fnname argnames arg indices
  deriving (Eq)
 
 instance Eq PolyExpr where
   (Plus []) == (Plus []) = True
   (Plus (x:xs)) == (Plus ys) =
-    case findIndex ((==) x) ys of
+    case elemIndex x ys of
       Just i -> let (hs, _:ts) = splitAt i ys in
-                  (Plus xs) == (Plus (hs ++ ts))
+                  Plus xs == Plus (hs ++ ts)
       Nothing -> False
   _ == _ = False
 
@@ -420,22 +430,22 @@
   (Term a []) == (Term b [])
     | a /= b =  False
     | otherwise = True
-  (Term a (((Quote x),n):xs)) == (Term b ys)
-    | (a /= b) && (a /= (negate b)) =  False
-    | otherwise = case findIndex ((==) ((Quote x),n)) ys of
+  (Term a ((Quote x, n):xs)) == (Term b ys)
+    | (a /= b) && (a /= negate b) =  False
+    | otherwise = case elemIndex (Quote x, n) ys of
                     Just i -> let (hs, _:ts) = splitAt i ys in
-                                (Term a xs) == (Term b (hs ++ ts))
-                    Nothing -> case findIndex ((==) ((Quote (mathNegate x)),n)) ys of
+                                Term a xs == Term b (hs ++ ts)
+                    Nothing -> case elemIndex (Quote (mathNegate x), n) ys of
                                  Just i -> let (hs, _:ts) = splitAt i ys in
                                              if even n
-                                               then (Term a xs) == (Term b (hs ++ ts))
-                                               else (Term (negate a) xs) == (Term b (hs ++ ts))
+                                               then Term a xs == Term b (hs ++ ts)
+                                               else Term (negate a) xs == Term b (hs ++ ts)
                                  Nothing -> False
   (Term a (x:xs)) == (Term b ys)
-    | (a /= b) && (a /= (negate b)) =  False
-    | otherwise = case findIndex ((==) x) ys of
+    | (a /= b) && (a /= negate b) =  False
+    | otherwise = case elemIndex x ys of
                     Just i -> let (hs, _:ts) = splitAt i ys in
-                                (Term a xs) == (Term b (hs ++ ts))
+                                Term a xs == Term b (hs ++ ts)
                     Nothing -> False
   _ == _ = False
 
@@ -449,7 +459,7 @@
   tensorElems :: a -> V.Vector a
   tensorSize :: a -> [Integer]
   tensorIndices :: a -> [Index EgisonValue]
-  fromTensor :: (Tensor a) -> EgisonM a
+  fromTensor :: Tensor a -> EgisonM a
   toTensor :: a -> EgisonM (Tensor a)
   undef :: a
 
@@ -457,7 +467,7 @@
   tensorElems (TensorData (Tensor _ xs _)) = xs
   tensorSize (TensorData (Tensor ns _ _)) = ns
   tensorIndices (TensorData (Tensor _ _ js)) = js
-  fromTensor t@(Tensor _ _ _) = return $ TensorData t
+  fromTensor t@Tensor{} = return $ TensorData t
   fromTensor (Scalar x) = return x
   toTensor (TensorData t) = return t
   toTensor x = return $ Scalar x
@@ -467,7 +477,7 @@
   tensorElems (Intermediate (ITensor (Tensor _ xs _))) = xs
   tensorSize (Intermediate (ITensor (Tensor ns _ _))) = ns
   tensorIndices (Intermediate (ITensor (Tensor _ _ js))) = js
-  fromTensor t@(Tensor _ _ _) = return $ Intermediate $ ITensor t
+  fromTensor t@Tensor{} = return $ Intermediate $ ITensor t
   fromTensor (Scalar x) = return x
   toTensor (Intermediate (ITensor t)) = return t
   toTensor x = return $ Scalar x
@@ -478,13 +488,16 @@
 --
 
 symbolScalarData :: String -> String -> EgisonValue
-symbolScalarData id name = ScalarData (Div (Plus [(Term 1 [(Symbol id name [], 1)])]) (Plus [(Term 1 [])]))
+symbolScalarData id name = ScalarData (Div (Plus [Term 1 [(Symbol id name [], 1)]]) (Plus [Term 1 []]))
 
 getSymId :: EgisonValue -> String
-getSymId (ScalarData (Div (Plus [(Term 1 [(Symbol id name [], 1)])]) (Plus [(Term 1 [])]))) = id
+getSymId (ScalarData (Div (Plus [Term 1 [(Symbol id name [], 1)]]) (Plus [Term 1 []]))) = id
 
+getSymName :: EgisonValue -> String
+getSymName (ScalarData (Div (Plus [Term 1 [(Symbol id name [], 1)]]) (Plus [Term 1 []]))) = name
+
 mathExprToEgison :: ScalarData -> EgisonValue
-mathExprToEgison (Div p1 p2) = InductiveData "Div" [(polyExprToEgison p1), (polyExprToEgison p2)]
+mathExprToEgison (Div p1 p2) = InductiveData "Div" [polyExprToEgison p1, polyExprToEgison p2]
 
 polyExprToEgison :: PolyExpr -> EgisonValue
 polyExprToEgison (Plus ts) = InductiveData "Plus" [Collection (Sq.fromList (map termExprToEgison ts))]
@@ -493,28 +506,43 @@
 termExprToEgison (Term a xs) = InductiveData "Term" [toEgison a, Collection (Sq.fromList (map symbolExprToEgison xs))]
 
 symbolExprToEgison :: (SymbolExpr, Integer) -> EgisonValue
-symbolExprToEgison (Symbol id x js, n) = Tuple [InductiveData "Symbol" [symbolScalarData id x, Collection (Sq.fromList (map (\j -> case j of
-                                                                                                                                     Superscript k -> InductiveData "Sup" [ScalarData k]
-                                                                                                                                     Subscript k -> InductiveData "Sub" [ScalarData k]
-                                                                                                                             ) js))], toEgison n]
+symbolExprToEgison (Symbol id x js, n) = Tuple [InductiveData "Symbol" [symbolScalarData id x, f js], toEgison n]
+ where
+  f js = Collection (Sq.fromList (map (\j -> case j of
+                                               Superscript k -> InductiveData "Sup" [ScalarData k]
+                                               Subscript k -> InductiveData "Sub" [ScalarData k]
+                                               Userscript k -> InductiveData "User" [ScalarData k]
+                                      ) js))
 symbolExprToEgison (Apply fn mExprs, n) = Tuple [InductiveData "Apply" [fn, Collection (Sq.fromList (map mathExprToEgison mExprs))], toEgison n]
 symbolExprToEgison (Quote mExpr, n) = Tuple [InductiveData "Quote" [mathExprToEgison mExpr], toEgison n]
+symbolExprToEgison (FunctionData name argnames args js, n) = case name of
+                                                               Nothing -> Tuple [InductiveData "Function" [symbolScalarData "" "", Collection (Sq.fromList argnames), Collection (Sq.fromList args), f js], toEgison n]
+                                                               Just name' -> Tuple [InductiveData "Function" [name', Collection (Sq.fromList argnames), Collection (Sq.fromList args), f js], toEgison n]
+ where
+  f js = Collection (Sq.fromList (map (\j -> case j of
+                                               Superscript k -> InductiveData "Sup" [ScalarData k]
+                                               Subscript k -> InductiveData "Sub" [ScalarData k]
+                                               Userscript k -> InductiveData "User" [ScalarData k]
+                                      ) js))
 
 egisonToScalarData :: EgisonValue -> EgisonM ScalarData
 egisonToScalarData (InductiveData "Div" [p1, p2]) = Div <$> egisonToPolyExpr p1 <*> egisonToPolyExpr p2
-egisonToScalarData p1@(InductiveData "Plus" _) = Div <$> egisonToPolyExpr p1 <*> (return (Plus [(Term 1 [])]))
+egisonToScalarData p1@(InductiveData "Plus" _) = Div <$> egisonToPolyExpr p1 <*> return (Plus [Term 1 []])
 egisonToScalarData t1@(InductiveData "Term" _) = do
   t1' <- egisonToTermExpr t1
-  return $ Div (Plus [t1']) (Plus [(Term 1 [])])
+  return $ Div (Plus [t1']) (Plus [Term 1 []])
 egisonToScalarData s1@(InductiveData "Symbol" _) = do
   s1' <- egisonToSymbolExpr (Tuple [s1, toEgison (1 ::Integer)])
-  return $ Div (Plus [(Term 1 [s1'])]) (Plus [(Term 1 [])])
+  return $ Div (Plus [Term 1 [s1']]) (Plus [Term 1 []])
 egisonToScalarData s1@(InductiveData "Apply" _) = do
   s1' <- egisonToSymbolExpr (Tuple [s1, toEgison (1 :: Integer)])
-  return $ Div (Plus [(Term 1 [s1'])]) (Plus [(Term 1 [])])
+  return $ Div (Plus [Term 1 [s1']]) (Plus [Term 1 []])
 egisonToScalarData s1@(InductiveData "Quote" _) = do
   s1' <- egisonToSymbolExpr (Tuple [s1, toEgison (1 :: Integer)])
-  return $ Div (Plus [(Term 1 [s1'])]) (Plus [(Term 1 [])])
+  return $ Div (Plus [Term 1 [s1']]) (Plus [Term 1 []])
+egisonToScalarData s1@(InductiveData "Function" _) = do
+  s1' <- egisonToSymbolExpr (Tuple [s1, toEgison (1 :: Integer)])
+  return $ Div (Plus [Term 1 [s1']]) (Plus [Term 1 []])
 egisonToScalarData val = liftError $ throwError $ TypeMismatch "math expression" (Value val)
 
 egisonToPolyExpr :: EgisonValue -> EgisonM PolyExpr
@@ -526,18 +554,19 @@
 egisonToTermExpr val = liftError $ throwError $ TypeMismatch "math term expression" (Value val)
 
 egisonToSymbolExpr :: EgisonValue -> EgisonM (SymbolExpr, Integer)
-egisonToSymbolExpr (Tuple [InductiveData "Symbol" [x, (Collection seq)], n]) = do
+egisonToSymbolExpr (Tuple [InductiveData "Symbol" [x, Collection seq], n]) = do
   let js = toList seq
   js' <- mapM (\j -> case j of
                          InductiveData "Sup" [ScalarData k] -> return (Superscript k)
                          InductiveData "Sub" [ScalarData k] -> return (Subscript k)
+                         InductiveData "User" [ScalarData k] -> return (Userscript k)
                          _ -> liftError $ throwError $ TypeMismatch "math symbol expression" (Value j)
                ) js
   n' <- fromEgison n
   case x of
-    (ScalarData (Div (Plus [(Term 1 [(Symbol id name [], 1)])]) (Plus [(Term 1 [])]))) ->
+    (ScalarData (Div (Plus [Term 1 [(Symbol id name [], 1)]]) (Plus [Term 1 []]))) ->
       return (Symbol id name js', n')
-egisonToSymbolExpr (Tuple [InductiveData "Apply" [fn, (Collection mExprs)], n]) = do
+egisonToSymbolExpr (Tuple [InductiveData "Apply" [fn, Collection mExprs], n]) = do
   mExprs' <- mapM egisonToScalarData (toList mExprs)
   n' <- fromEgison n
   return (Apply fn mExprs', n')
@@ -545,6 +574,19 @@
   mExpr' <- egisonToScalarData mExpr
   n' <- fromEgison n
   return (Quote mExpr', n')
+egisonToSymbolExpr (Tuple [InductiveData "Function" [name, (Collection argnames), (Collection args), Collection seq], n]) = do
+  let js = toList seq
+  js' <- mapM (\j -> case j of
+                         InductiveData "Sup" [ScalarData k] -> return (Superscript k)
+                         InductiveData "Sub" [ScalarData k] -> return (Subscript k)
+                         InductiveData "User" [ScalarData k] -> return (Userscript k)
+                         _ -> liftError $ throwError $ TypeMismatch "math symbol expression" (Value j)
+               ) js
+  n' <- fromEgison n
+  let name' = case getSymName name of
+                "" -> Nothing
+                s -> Just $ name
+  return (FunctionData name' (toList argnames) (toList args) js', n')
 egisonToSymbolExpr val = liftError $ throwError $ TypeMismatch "math symbol expression" (Value val)
 
 mathNormalize' :: ScalarData -> ScalarData
@@ -555,32 +597,31 @@
  where
   f :: TermExpr -> [TermExpr] -> TermExpr
   f ret [] =  ret
-  f (Term a xs) ((Term b ys):ts) =
+  f (Term a xs) (Term b ys:ts) =
     f (Term (gcd a b) (g xs ys)) ts
   g :: [(SymbolExpr, Integer)] -> [(SymbolExpr, Integer)] -> [(SymbolExpr, Integer)]
   g [] ys = []
   g ((x, n):xs) ys = let (z, m) = h (x, n) ys in
-    if m == 0 then g xs ys else ((z, m):(g xs ys))
+    if m == 0 then g xs ys else (z, m):g xs ys
   h :: (SymbolExpr, Integer) -> [(SymbolExpr, Integer)] -> (SymbolExpr, Integer)
   h (x, n) [] = (x, 0)
-  h ((Quote x), n) (((Quote y), m):ys) = if x == y
-                                        then ((Quote x), (min n m))
-                                        else if x == (mathNegate y)
-                                             then ((Quote x), (min n m))
-                                             else h ((Quote x), n) ys
+  h (Quote x, n) ((Quote y, m):ys)
+    | x == y = (Quote x, min n m)
+    | x == mathNegate y = (Quote x, min n m)
+    | otherwise = h (Quote x, n) ys
   h (x, n) ((y, m):ys) = if x == y
-                         then (x, (min n m))
+                         then (x, min n m)
                          else h (x, n) ys
 
 mathDivide :: ScalarData -> ScalarData
-mathDivide (Div (Plus ts1) (Plus [])) = (Div (Plus ts1) (Plus []))
-mathDivide (Div (Plus []) (Plus ts2)) = (Div (Plus []) (Plus ts2))
+mathDivide (Div (Plus ts1) (Plus [])) = Div (Plus ts1) (Plus [])
+mathDivide (Div (Plus []) (Plus ts2)) = Div (Plus []) (Plus ts2)
 mathDivide (Div (Plus ts1) (Plus ts2)) =
   let z = termsGcd (ts1 ++ ts2) in
   case z of
 --    (Term 1 []) -> (Div (Plus ts1) (Plus ts2))
     (Term c zs) -> case ts2 of
-      [(Term a _)] -> if a < 0
+      [Term a _] -> if a < 0
                       then (Div (Plus (map (\t -> mathDivideTerm t (Term (-1 * c) zs)) ts1)) (Plus (map (\t -> mathDivideTerm t (Term (-1 * c) zs)) ts2)))
                       else (Div (Plus (map (\t -> mathDivideTerm t z) ts1)) (Plus (map (\t -> mathDivideTerm t z) ts2)))
       _ -> (Div (Plus (map (\t -> mathDivideTerm t z) ts1)) (Plus (map (\t -> mathDivideTerm t z) ts2)))
@@ -588,23 +629,21 @@
 mathDivideTerm :: TermExpr -> TermExpr -> TermExpr
 mathDivideTerm (Term a xs) (Term b ys) =
   let (sgn, zs) = f 1 xs ys in
-  (Term (sgn * (div a b)) zs)
+  Term (sgn * div a b) zs
  where
   f :: Integer -> [(SymbolExpr, Integer)] -> [(SymbolExpr, Integer)] -> (Integer, [(SymbolExpr, Integer)])
   f sgn xs [] = (sgn, xs)
   f sgn xs ((y, n):ys) =
     let (sgns, zs) = unzip (map (\(x, m) -> (g (x, m) (y, n))) xs) in
-    f (sgn * (product sgns)) zs ys
+    f (sgn * product sgns) zs ys
   g :: (SymbolExpr, Integer) -> (SymbolExpr, Integer) -> (Integer, (SymbolExpr, Integer))
-  g ((Quote x), n) ((Quote y), m) =
-    if x == y
-    then (1, ((Quote x), (n - m)))
-    else if x == (mathNegate y)
-         then if even m then (1, ((Quote x), (n - m))) else (-1, ((Quote x), (n - m)))
-         else (1, ((Quote x), n))
+  g (Quote x, n) (Quote y, m)
+    | x == y = (1, (Quote x, n - m))
+    | x == mathNegate y = if even m then (1, (Quote x, n - m)) else (-1, (Quote x, n - m))
+    | otherwise = (1, (Quote x, n))
   g (x, n) (y, m) =
     if x == y
-    then (1, (x, (n - m)))
+    then (1, (x, n - m))
     else (1, (x, n))
               
 mathRemoveZeroSymbol :: ScalarData -> ScalarData
@@ -625,29 +664,28 @@
       _ -> Div (Plus ts1') (Plus ts2')
 
 mathFold :: ScalarData -> ScalarData
-mathFold mExpr = (mathTermFold (mathSymbolFold (mathTermFold mExpr)))
+mathFold mExpr = mathTermFold (mathSymbolFold (mathTermFold mExpr))
 
 mathSymbolFold :: ScalarData -> ScalarData
 mathSymbolFold (Div (Plus ts1) (Plus ts2)) = Div (Plus (map f ts1)) (Plus (map f ts2))
  where
   f :: TermExpr -> TermExpr
   f (Term a xs) = let (ys, sgns) = unzip $ g [] xs
-                    in Term ((product sgns) * a) ys
+                    in Term (product sgns * a) ys
   g :: [((SymbolExpr, Integer),Integer)] -> [(SymbolExpr, Integer)] -> [((SymbolExpr, Integer),Integer)]
   g ret [] = ret
   g ret ((x, n):xs) =
-    if (any (p (x, n)) ret)
+    if any (p (x, n)) ret
       then g (map (h (x, n)) ret) xs
       else g (ret ++ [((x, n), 1)]) xs
   p :: (SymbolExpr, Integer) -> ((SymbolExpr, Integer), Integer) -> Bool
-  p ((Quote x), _) (((Quote y), _),_) = (x == y) || ((mathNegate x) == y)
+  p (Quote x, _) ((Quote y, _),_) = (x == y) || (mathNegate x == y)
   p (x, _) ((y, _),_) = x == y
   h :: (SymbolExpr, Integer) -> ((SymbolExpr, Integer), Integer) -> ((SymbolExpr, Integer), Integer)
-  h ((Quote x), n) (((Quote y), m), sgn) = if x == y
-                                      then (((Quote y), m + n), sgn)
-                                      else if x == (mathNegate y)
-                                            then if even n then (((Quote y), m + n), sgn) else (((Quote y), m + n), -1 * sgn)
-                                            else (((Quote y), m), sgn)
+  h (Quote x, n) ((Quote y, m), sgn)
+    | x == y = ((Quote y, m + n), sgn)
+    | x == mathNegate y = if even n then ((Quote y, m + n), sgn) else ((Quote y, m + n), -1 * sgn)
+    | otherwise = ((Quote y, m), sgn)
   h (x, n) ((y, m), sgn) = if x == y
                              then ((y, m + n), sgn)
                              else ((y, m), sgn)
@@ -656,17 +694,17 @@
 mathTermFold (Div (Plus ts1) (Plus ts2)) = Div (Plus (f ts1)) (Plus (f ts2))
  where
   f :: [TermExpr] -> [TermExpr]
-  f ts = f' [] ts
+  f = f' []
   f' :: [TermExpr] -> [TermExpr] -> [TermExpr]
   f' ret [] = ret
-  f' ret ((Term a xs):ts) =
+  f' ret (Term a xs:ts) =
     if any (\(Term _ ys) -> (fst (p 1 xs ys))) ret
       then f' (map (g (Term a xs)) ret) ts
-      else f' (ret ++ [(Term a xs)]) ts
+      else f' (ret ++ [Term a xs]) ts
   g :: TermExpr -> TermExpr -> TermExpr
   g (Term a xs) (Term b ys) = let (c, sgn) = p 1 xs ys in
                                 if c
-                                  then (Term ((sgn * a) + b) ys)
+                                  then Term ((sgn * a) + b) ys
                                   else Term b ys
   p :: Integer -> [(SymbolExpr, Integer)] -> [(SymbolExpr, Integer)] -> (Bool, Integer)
   p sgn [] [] = (True, sgn)
@@ -678,14 +716,13 @@
         else (False, 0)
   q :: (SymbolExpr, Integer) -> [(SymbolExpr, Integer)] -> [(SymbolExpr, Integer)] -> (Bool, [(SymbolExpr, Integer)], Integer)
   q _ _ [] = (False, [], 1)
-  q ((Quote x), n) ret (((Quote y), m):ys) = if (x == y) && (n == m)
-                                               then (True, (ret ++ ys), 1)
-                                               else if ((mathNegate x) == y) && (n == m)
-                                                      then if even n then (True, (ret ++ ys), 1) else (True, (ret ++ ys), -1)
-                                                      else q ((Quote x), n) (ret ++ [((Quote y), m)]) ys
-  q ((Quote x), n) ret ((y,m):ys) = q ((Quote x), n) (ret ++ [(y, m)]) ys
+  q (Quote x, n) ret ((Quote y, m):ys)
+    | (x == y) && (n == m) = (True, ret ++ ys, 1)
+    | (mathNegate x == y) && (n == m) = if even n then (True, ret ++ ys, 1) else (True, ret ++ ys, -1)
+    | otherwise = q (Quote x, n) (ret ++ [(Quote y, m)]) ys
+  q (Quote x, n) ret ((y,m):ys) = q (Quote x, n) (ret ++ [(y, m)]) ys
   q (x, n) ret ((y, m):ys) = if (x == y) && (n == m)
-                               then (True, (ret ++ ys), 1)
+                               then (True, ret ++ ys, 1)
                                else q (x, n) (ret ++ [(y, m)]) ys
 
 --
@@ -716,10 +753,10 @@
 mathNegate' (Plus ts) = Plus (map (\(Term a xs) -> (Term (negate a) xs)) ts)
 
 mathNumerator :: ScalarData -> ScalarData
-mathNumerator (Div m _) = Div m (Plus [(Term 1 [])])
+mathNumerator (Div m _) = Div m (Plus [Term 1 []])
 
 mathDenominator :: ScalarData -> ScalarData
-mathDenominator (Div _ n) = Div n (Plus [(Term 1 [])])
+mathDenominator (Div _ n) = Div n (Plus [Term 1 []])
 
 --
 -- ExtractScalar
@@ -731,33 +768,35 @@
 
 extractScalar' :: WHNFData -> EgisonM ScalarData
 extractScalar' (Value (ScalarData x)) = return x
-extractScalar' val = throwError $ TypeMismatch "integer or string" $ val
+extractScalar' val = throwError $ TypeMismatch "integer or string" val
 
 --
 -- Tensors
 --
 
-initTensor :: [Integer] -> [a] -> [EgisonValue] -> [EgisonValue] -> (Tensor a)
-initTensor ns xs sup sub = Tensor ns (V.fromList xs) ((map Superscript sup) ++ (map Subscript sub))
+initTensor :: [Integer] -> [a] -> [EgisonValue] -> [EgisonValue] -> Tensor a
+initTensor ns xs sup sub = Tensor ns (V.fromList xs) (map Superscript sup ++ map Subscript sub)
 
-tSize :: (Tensor a) -> [Integer]
+tSize :: Tensor a -> [Integer]
 tSize (Tensor ns _ _) = ns
 tSize (Scalar _) = []
 
-tToList :: (Tensor a) -> [a]
+tToList :: Tensor a -> [a]
 tToList (Tensor _ xs _) = V.toList xs
 tToList (Scalar x) = [x]
 
-tToVector :: (Tensor a) -> V.Vector a
+tToVector :: Tensor a -> V.Vector a
 tToVector (Tensor _ xs _) = xs
 tToVector (Scalar x) = V.fromList [x]
 
-tIndex :: (Tensor a) -> [Index EgisonValue]
+tIndex :: Tensor a -> [Index EgisonValue]
 tIndex (Tensor _ _ js) = js
 tIndex (Scalar _) = []
 
-tIntRef' :: HasTensor a => Integer -> (Tensor a) -> EgisonM a
-tIntRef' i (Tensor [_] xs _) = fromTensor $ Scalar $ xs V.! (fromIntegral (i - 1))
+tIntRef' :: HasTensor a => Integer -> Tensor a -> EgisonM a
+tIntRef' i (Tensor [ary] xs _) = let n = fromIntegral (length [ary]) in
+                                     if (0 < i) && (i <= (n + 3)) then fromTensor $ Scalar $ xs V.! fromIntegral (i - 1)
+                                                                  else throwError $ TensorIndexOutOfBounds i (n + 3)
 tIntRef' i (Tensor (n:ns) xs js) =
   if (0 < i) && (i <= n)
    then let w = fromIntegral (product ns) in
@@ -766,47 +805,50 @@
    else throwError $ TensorIndexOutOfBounds i n
 tIntRef' i _ = throwError $ Default "More indices than the order of the tensor"
  
-tIntRef :: HasTensor a => [Integer] -> (Tensor a) -> EgisonM (Tensor a)
+tIntRef :: HasTensor a => [Integer] -> Tensor a -> EgisonM (Tensor a)
 tIntRef [] (Tensor [] xs _)
   | V.length xs == 1 = return $ Scalar (xs V.! 0)
   | otherwise = throwError $ EgisonBug "sevaral elements in scalar tensor"
 tIntRef [] t = return t
 tIntRef (m:ms) t = tIntRef' m t >>= toTensor >>= tIntRef ms 
 
-tref :: HasTensor a => [Index EgisonValue] -> (Tensor a) -> EgisonM a
+tref :: HasTensor a => [Index EgisonValue] -> Tensor a -> EgisonM a
 tref [] (Tensor [] xs _)
   | V.length xs == 1 = fromTensor $ Scalar (xs V.! 0)
   | otherwise = throwError $ EgisonBug "sevaral elements in scalar tensor"
 tref [] t = fromTensor t
-tref ((Subscript (ScalarData (Div (Plus [(Term m [])]) (Plus [(Term 1 [])])))):ms) t = tIntRef' m t >>= toTensor >>= tref ms
-tref ((Superscript (ScalarData (Div (Plus [(Term m [])]) (Plus [(Term 1 [])])))):ms) t = tIntRef' m t >>= toTensor >>= tref ms
-tref ((SupSubscript (ScalarData (Div (Plus [(Term m [])]) (Plus [(Term 1 [])])))):ms) t = tIntRef' m t >>= toTensor >>= tref ms
-tref ((Subscript (Tuple [mVal, nVal])):ms) t@(Tensor is _ _) = do
+tref (Subscript (ScalarData (Div (Plus [Term m []]) (Plus [Term 1 []]))):ms) t = tIntRef' m t >>= toTensor >>= tref ms
+tref (Subscript (ScalarData (Div (Plus []) (Plus [Term 1 []]))):ms) t = tIntRef' 0 t >>= toTensor >>= tref ms
+tref (Superscript (ScalarData (Div (Plus [Term m []]) (Plus [Term 1 []]))):ms) t = tIntRef' m t >>= toTensor >>= tref ms
+tref (Superscript (ScalarData (Div (Plus []) (Plus [Term 1 []]))):ms) t = tIntRef' 0 t >>= toTensor >>= tref ms
+tref (SupSubscript (ScalarData (Div (Plus [Term m []]) (Plus [Term 1 []]))):ms) t = tIntRef' m t >>= toTensor >>= tref ms
+tref (SupSubscript (ScalarData (Div (Plus []) (Plus [Term 1 []]))):ms) t = tIntRef' 0 t >>= toTensor >>= tref ms
+tref (Subscript (Tuple [mVal, nVal]):ms) t@(Tensor is _ _) = do
   m <- fromEgison mVal
   n <- fromEgison nVal
   if m > n
-    then do
-      fromTensor (Tensor (take (length is) (repeat 0)) V.empty [])
+    then
+      fromTensor (Tensor (replicate (length is) 0) V.empty [])
     else do
       ts <- mapM (\i -> tIntRef' i t >>= toTensor >>= tref ms >>= toTensor) [m..n]
       symId <- fresh
       tConcat (Subscript (symbolScalarData "" (":::" ++ symId))) ts >>= fromTensor
-tref ((Superscript (Tuple [mVal, nVal])):ms) t@(Tensor is _ _) = do
+tref (Superscript (Tuple [mVal, nVal]):ms) t@(Tensor is _ _) = do
   m <- fromEgison mVal
   n <- fromEgison nVal
   if m > n
-    then do
-      fromTensor (Tensor (take (length is) (repeat 0)) V.empty [])
+    then
+      fromTensor (Tensor (replicate (length is) 0) V.empty [])
     else do
       ts <- mapM (\i -> tIntRef' i t >>= toTensor >>= tref ms >>= toTensor) [m..n]
       symId <- fresh
       tConcat (Superscript (symbolScalarData "" (":::" ++ symId))) ts >>= fromTensor
-tref ((SupSubscript (Tuple [mVal, nVal])):ms) t@(Tensor is _ _) = do
+tref (SupSubscript (Tuple [mVal, nVal]):ms) t@(Tensor is _ _) = do
   m <- fromEgison mVal
   n <- fromEgison nVal
   if m > n
-    then do
-      fromTensor (Tensor (take (length is) (repeat 0)) V.empty [])
+    then
+      fromTensor (Tensor (replicate (length is) 0) V.empty [])
     else do
       ts <- mapM (\i -> tIntRef' i t >>= toTensor >>= tref ms >>= toTensor) [m..n]
       symId <- fresh
@@ -819,44 +861,43 @@
 
 enumTensorIndices :: [Integer] -> [[Integer]]
 enumTensorIndices [] = [[]]
-enumTensorIndices (n:ns) = concat (map (\i -> (map (\is -> i:is) (enumTensorIndices ns))) [1..n])
+enumTensorIndices (n:ns) = concatMap (\i -> (map (\is -> i:is) (enumTensorIndices ns))) [1..n]
 
 transIndex :: [Index EgisonValue] -> [Index EgisonValue] -> [Integer] -> EgisonM [Integer]
 transIndex [] [] is = return is
 transIndex (j1:js1) js2 is = do
   let (hjs2, tjs2) = break (\j2 -> j1 == j2) js2
   if tjs2 == []
-    then do throwError $ InconsistentTensorIndex
-    else do let n = (length hjs2) + 1
-            rs <- transIndex js1 (hjs2 ++ (tail tjs2)) ((take (n - 1) is) ++ (drop n is))
-            return ((nth (fromIntegral n) is):rs)
-transIndex _ _ _ = throwError $ InconsistentTensorSize
+    then throwError InconsistentTensorIndex
+    else do let n = length hjs2 + 1
+            rs <- transIndex js1 (hjs2 ++ tail tjs2) (take (n - 1) is ++ drop n is)
+            return (nth (fromIntegral n) is:rs)
+transIndex _ _ _ = throwError InconsistentTensorSize
 
-tTranspose :: HasTensor a => [Index EgisonValue] -> (Tensor a) -> EgisonM (Tensor a)
+tTranspose :: HasTensor a => [Index EgisonValue] -> Tensor a -> EgisonM (Tensor a)
 tTranspose is t@(Tensor ns xs js) = do
   ns' <- transIndex js is ns
-  xs' <- mapM (transIndex js is) (enumTensorIndices ns') >>= mapM (flip tIntRef t) >>= mapM fromTensor >>= return . V.fromList
+  xs' <- mapM (transIndex js is) (enumTensorIndices ns') >>= mapM (`tIntRef` t) >>= mapM fromTensor >>= return . V.fromList
   return $ Tensor ns' xs' is
 
-tTranspose' :: HasTensor a => [EgisonValue] -> (Tensor a) -> EgisonM (Tensor a)
+tTranspose' :: HasTensor a => [EgisonValue] -> Tensor a -> EgisonM (Tensor a)
 tTranspose' is t@(Tensor ns xs js) = do
   is' <- g is js
   tTranspose is' t
  where
-  f :: (Index EgisonValue) -> EgisonValue
+  f :: Index EgisonValue -> EgisonValue
   f (Subscript i) = i
   f (Superscript i) = i
   f (SupSubscript i) = i
   g :: [EgisonValue] -> [Index EgisonValue] -> EgisonM [Index EgisonValue]
   g [] js = return []
-  g (i:is) js = case find (\j -> i == (f j)) js of
-                  Nothing ->  throwError $ InconsistentTensorIndex
+  g (i:is) js = case find (\j -> i == f j) js of
+                  Nothing ->  throwError InconsistentTensorIndex
                   (Just j') -> do js' <- g is js
                                   return $ j':js'
 
-tFlipIndices :: HasTensor a => (Tensor a) -> EgisonM (Tensor a)
-tFlipIndices (Tensor ns xs js) = do
-  return $ Tensor ns xs (map flipIndex js)
+tFlipIndices :: HasTensor a => Tensor a -> EgisonM (Tensor a)
+tFlipIndices (Tensor ns xs js) = return $ Tensor ns xs (map flipIndex js)
  where
   flipIndex (Subscript i) = Superscript i
   flipIndex (Superscript i) = Subscript i
@@ -864,11 +905,11 @@
 
 appendDFscripts :: Integer -> WHNFData -> EgisonM WHNFData
 appendDFscripts id (Intermediate (ITensor (Tensor s xs is))) = do
-  let k = fromIntegral ((length s) - (length is))
-  return $ Intermediate (ITensor (Tensor s xs (is ++ (map (DFscript id) [1..k]))))
+  let k = fromIntegral (length s - length is)
+  return $ Intermediate (ITensor (Tensor s xs (is ++ map (DFscript id) [1..k])))
 appendDFscripts id (Value (TensorData (Tensor s xs is))) = do
-  let k = fromIntegral ((length s) - (length is))
-  return $ Value (TensorData (Tensor s xs (is ++ (map (DFscript id) [1..k]))))
+  let k = fromIntegral (length s - length is)
+  return $ Value (TensorData (Tensor s xs (is ++ map (DFscript id) [1..k])))
 appendDFscripts _ whnf = return whnf
 
 removeDFscripts :: WHNFData -> EgisonM WHNFData
@@ -888,108 +929,108 @@
   isDF _ = False
 removeDFscripts whnf = return whnf
 
-tMap :: HasTensor a => (a -> EgisonM a) -> (Tensor a) -> EgisonM (Tensor a)
+tMap :: HasTensor a => (a -> EgisonM a) -> Tensor a -> EgisonM (Tensor a)
 tMap f (Tensor ns xs js') = do
-  let k = fromIntegral $ (length ns) - (length js')
-  let js = (js' ++ (map (DFscript 0) [1..k]))
+  let k = fromIntegral $ length ns - length js'
+  let js = js' ++ map (DFscript 0) [1..k]
   xs' <- mapM f (V.toList xs) >>= return . V.fromList
   t <- toTensor (V.head xs')
   case t of
     (Tensor ns1 _ js1') -> do
-      let k1 = fromIntegral $ (length ns1) - (length js1')
-      let js1 = (js1' ++ (map (DFscript 0) [1..k1]))
+      let k1 = fromIntegral $ length ns1 - length js1'
+      let js1 = js1' ++ map (DFscript 0) [1..k1]
       tContract' $ Tensor (ns ++ ns1) (V.concat (V.toList (V.map tensorElems xs'))) (js ++ js1)
     _ -> return $ Tensor ns xs' js
-tMap f (Scalar x) = f x >>= return . Scalar
+tMap f (Scalar x) = Scalar <$> f x
 
 tMapN :: HasTensor a => ([a] -> EgisonM a) -> [Tensor a] -> EgisonM (Tensor a)
-tMapN f ts@((Tensor ns xs js):_) = do
+tMapN f ts@(Tensor ns xs js:_) = do
   xs' <- mapM (\is -> mapM (tIntRef is) ts >>= mapM fromTensor >>= f) (enumTensorIndices ns)
   return $ Tensor ns (V.fromList xs') js
-tMapN f xs = mapM fromTensor xs >>= f >>= return . Scalar
+tMapN f xs = Scalar <$> (mapM fromTensor xs >>= f)
 
 tMap2 :: HasTensor a => (a -> a -> EgisonM a) -> Tensor a -> Tensor a -> EgisonM (Tensor a)
 tMap2 f t1@(Tensor ns1 xs1 js1') t2@(Tensor ns2 xs2 js2') = do
-  let k1 = fromIntegral $ (length ns1) - (length js1')
-  let js1 = (js1' ++ (map (DFscript 0) [1..k1]))
-  let k2 = fromIntegral $ (length ns2) - (length js2')
-  let js2 = (js2' ++ (map (DFscript 0) [1..k2]))
+  let k1 = fromIntegral $ length ns1 - length js1'
+  let js1 = js1' ++ map (DFscript 0) [1..k1]
+  let k2 = fromIntegral $ length ns2 - length js2'
+  let js2 = js2' ++ map (DFscript 0) [1..k2]
   let (cjs, tjs1, tjs2) = h js1 js2
   t1' <- tTranspose (cjs ++ tjs1) t1
   t2' <- tTranspose (cjs ++ tjs2) t2
   let cns = take (length cjs) (tSize t1')
-  rts1 <- mapM (flip tIntRef t1') (enumTensorIndices cns)
-  rts2 <- mapM (flip tIntRef t2') (enumTensorIndices cns)
+  rts1 <- mapM (`tIntRef` t1') (enumTensorIndices cns)
+  rts2 <- mapM (`tIntRef` t2') (enumTensorIndices cns)
   rts' <- mapM (\(t1, t2) -> tProduct f t1 t2) (zip rts1 rts2)
-  let ret = Tensor (cns ++ (tSize (head rts'))) (V.concat (map tToVector rts')) (cjs ++ tIndex (head rts'))
+  let ret = Tensor (cns ++ tSize (head rts')) (V.concat (map tToVector rts')) (cjs ++ tIndex (head rts'))
   tTranspose (uniq (tDiagIndex (js1 ++ js2))) ret
  where
   h :: [Index EgisonValue] -> [Index EgisonValue] -> ([Index EgisonValue], [Index EgisonValue], [Index EgisonValue])
-  h js1 js2 = let cjs = filter (\j -> elem j js2) js1 in
+  h js1 js2 = let cjs = filter (`elem` js2) js1 in
                 (cjs, js1 \\ cjs, js2 \\ cjs)
   uniq :: [Index EgisonValue] -> [Index EgisonValue]
   uniq [] = []
-  uniq (x:xs) = x:(uniq (delete x xs))
-tMap2 f t@(Tensor _ _ _) (Scalar x) = tMap (flip f x) t
+  uniq (x:xs) = x:uniq (delete x xs)
+tMap2 f t@(Tensor _ _ _) (Scalar x) = tMap (`f` x) t
 tMap2 f (Scalar x) t@(Tensor _ _ _) = tMap (f x) t
-tMap2 f (Scalar x1) (Scalar x2) = f x1 x2 >>= return . Scalar
+tMap2 f (Scalar x1) (Scalar x2) = Scalar <$> f x1 x2
 
 tDiag :: HasTensor a => Tensor a -> EgisonM (Tensor a)
 tDiag t@(Tensor _ _ js) = do
   case filter (\j -> any (p j) js) js of
     [] -> return t
     xs -> do
-      let ys = js \\ (xs ++ (map rev xs))
-      t2 <- tTranspose (xs ++ (map rev xs) ++ ys) t
+      let ys = js \\ (xs ++ map rev xs)
+      t2 <- tTranspose (xs ++ map rev xs ++ ys) t
       let (ns1, tmp) = splitAt (length xs) (tSize t2)
       let (_, ns2) = splitAt (length xs) tmp
       ts <- mapM (\is -> tIntRef (is ++ is) t2) (enumTensorIndices ns1)
-      return $ Tensor (ns1 ++ ns2) (V.concat (map tToVector ts)) ((map g xs) ++ ys)
+      return $ Tensor (ns1 ++ ns2) (V.concat (map tToVector ts)) (map g xs ++ ys)
  where
   p :: Index EgisonValue -> Index EgisonValue -> Bool
   p (Superscript i) (Subscript j) = i == j
   p (Subscript i) _ = False
   p _ _ = False
   rev :: Index EgisonValue -> Index EgisonValue
-  rev (Superscript i) = (Subscript i)
-  rev (Subscript i) = (Superscript i)
+  rev (Superscript i) = Subscript i
+  rev (Subscript i) = Superscript i
   g :: Index EgisonValue -> Index EgisonValue
-  g (Superscript i) = (SupSubscript i)
-  g (Subscript i) = (SupSubscript i)
+  g (Superscript i) = SupSubscript i
+  g (Subscript i) = SupSubscript i
 tDiag t = return t
 
 tDiagIndex :: [Index EgisonValue] -> [Index EgisonValue]
 tDiagIndex js =
   let xs = filter (\j -> any (p j) js) js in
-  let ys = js \\ (xs ++ (map rev xs)) in
-    (map g xs) ++ ys
+  let ys = js \\ (xs ++ map rev xs) in
+    map g xs ++ ys
  where
   p :: Index EgisonValue -> Index EgisonValue -> Bool
   p (Superscript i) (Subscript j) = i == j
   p (Subscript _) _ = False
   p _ _ = False
   rev :: Index EgisonValue -> Index EgisonValue
-  rev (Superscript i) = (Subscript i)
-  rev (Subscript i) = (Superscript i)
+  rev (Superscript i) = Subscript i
+  rev (Subscript i) = Superscript i
   g :: Index EgisonValue -> Index EgisonValue
-  g (Superscript i) = (SupSubscript i)
-  g (Subscript i) = (SupSubscript i)
+  g (Superscript i) = SupSubscript i
+  g (Subscript i) = SupSubscript i
 
-tSum :: HasTensor a => (a -> a -> EgisonM a) -> (Tensor a) -> (Tensor a) -> EgisonM (Tensor a)
+tSum :: HasTensor a => (a -> a -> EgisonM a) -> Tensor a -> Tensor a -> EgisonM (Tensor a)
 tSum f t1@(Tensor ns1 xs1 js1) t2@(Tensor _ _ _) = do
   t2' <- tTranspose js1 t2
   case t2' of
     (Tensor ns2 xs2 _)
       | ns2 == ns1 -> do ys <- V.mapM (\(x1,x2) -> f x1 x2) (V.zip xs1 xs2)
                          return (Tensor ns1 ys js1)
-      | otherwise -> throwError $ InconsistentTensorSize
+      | otherwise -> throwError InconsistentTensorSize
 
-tProduct :: HasTensor a => (a -> a -> EgisonM a) -> (Tensor a) -> (Tensor a) -> EgisonM (Tensor a)
+tProduct :: HasTensor a => (a -> a -> EgisonM a) -> Tensor a -> Tensor a -> EgisonM (Tensor a)
 tProduct f t1''@(Tensor ns1 xs1 js1') t2''@(Tensor ns2 xs2 js2') = do
-  let k1 = fromIntegral $ (length ns1) - (length js1')
-  let js1 = (js1' ++ (map (DFscript 0) [1..k1]))
-  let k2 = fromIntegral $ (length ns2) - (length js2')
-  let js2 = (js2' ++ (map (DFscript 0) [1..k2]))
+  let k1 = fromIntegral $ length ns1 - length js1'
+  let js1 = js1' ++ map (DFscript 0) [1..k1]
+  let k2 = fromIntegral $ length ns2 - length js2'
+  let js2 = js2' ++ map (DFscript 0) [1..k2]
   let (cjs1, cjs2, tjs1, tjs2) = h js1 js2
   let t1 = (Tensor ns1 xs1 js1)
   let t2 = (Tensor ns2 xs2 js2)
@@ -1009,45 +1050,45 @@
       rts' <- mapM (\is -> do rt1 <- tIntRef is t1'
                               rt2 <- tIntRef is t2'
                               tProduct f rt1 rt2) (enumTensorIndices cns1)
-      let ret = Tensor (cns1 ++ (tSize (head rts'))) (V.concat (map tToVector rts')) ((map g cjs1) ++ tIndex (head rts'))
-      ret2 <- tTranspose (uniq ((map g cjs1) ++ tjs1 ++ tjs2)) ret
+      let ret = Tensor (cns1 ++ tSize (head rts')) (V.concat (map tToVector rts')) (map g cjs1 ++ tIndex (head rts'))
+      ret2 <- tTranspose (uniq (map g cjs1 ++ tjs1 ++ tjs2)) ret
       return ret2
  where
   h :: [Index EgisonValue] -> [Index EgisonValue] -> ([Index EgisonValue], [Index EgisonValue], [Index EgisonValue], [Index EgisonValue])
   h js1 js2 = let cjs = filter (\j -> any (p j) js2) js1 in
-                (cjs, map rev cjs, js1 \\ cjs, js2 \\ (map rev cjs))
+                (cjs, map rev cjs, js1 \\ cjs, js2 \\ map rev cjs)
   p :: Index EgisonValue -> Index EgisonValue -> Bool
   p (Superscript i) (Subscript j) = i == j
   p (Subscript i) (Superscript j) = i == j
   p _ _ = False
   rev :: Index EgisonValue -> Index EgisonValue
-  rev (Superscript i) = (Subscript i)
-  rev (Subscript i) = (Superscript i)
+  rev (Superscript i) = Subscript i
+  rev (Subscript i) = Superscript i
   g :: Index EgisonValue -> Index EgisonValue
-  g (Superscript i) = (SupSubscript i)
-  g (Subscript i) = (SupSubscript i)
+  g (Superscript i) = SupSubscript i
+  g (Subscript i) = SupSubscript i
   uniq :: [Index EgisonValue] -> [Index EgisonValue]
   uniq [] = []
-  uniq (x:xs) = x:(uniq (delete x xs))
+  uniq (x:xs) = x:uniq (delete x xs)
 tProduct f (Scalar x) (Tensor ns xs js) = do
   xs' <- V.mapM (f x) xs
   return $ Tensor ns xs' js
 tProduct f (Tensor ns xs js) (Scalar x) = do
-  xs' <- V.mapM (flip f x) xs
+  xs' <- V.mapM (`f` x) xs
   return $ Tensor ns xs' js
-tProduct f (Scalar x1) (Scalar x2) = f x1 x2 >>= return . Scalar
+tProduct f (Scalar x1) (Scalar x2) = Scalar <$> f x1 x2
 
-tContract :: HasTensor a => (Tensor a) -> EgisonM [Tensor a]
+tContract :: HasTensor a => Tensor a -> EgisonM [Tensor a]
 tContract t = do
   t' <- tDiag t
   case t' of
-    (Tensor (n:ns) xs ((SupSubscript i):js)) -> do
-      ts <- mapM (\k -> tIntRef' k t') [1..n]
+    (Tensor (n:ns) xs (SupSubscript i:js)) -> do
+      ts <- mapM (`tIntRef'` t') [1..n]
       tss <- mapM toTensor ts >>= mapM tContract
       return $ concat tss
     _ -> return [t']
 
-tContract' :: HasTensor a => (Tensor a) -> EgisonM (Tensor a)
+tContract' :: HasTensor a => Tensor a -> EgisonM (Tensor a)
 tContract' t@(Tensor ns xs js) = do
   case findPairs p js of
     [] -> return t
@@ -1055,8 +1096,8 @@
       let ns' = (ns !! m):removePairs (m,n) ns
       let js' = (js !! m):removePairs (m,n) js
       let (hjs, mjs, tjs) = removePairs' (m,n) js
-      xs' <- mapM (\i -> (tref (hjs ++ [Subscript (ScalarData (Div (Plus [(Term i [])]) (Plus [(Term 1 [])])))] ++ mjs
-                                    ++ [Subscript (ScalarData (Div (Plus [(Term i [])]) (Plus [(Term 1 [])])))] ++ tjs) t))
+      xs' <- mapM (\i -> (tref (hjs ++ [Subscript (ScalarData (Div (Plus [Term i []]) (Plus [Term 1 []])))] ++ mjs
+                                    ++ [Subscript (ScalarData (Div (Plus [Term i []]) (Plus [Term 1 []])))] ++ tjs) t))
                   [1..(ns !! m)]
       mapM toTensor xs' >>= tConcat (js !! m) >>= tTranspose (hjs ++ [js !! m] ++ mjs ++ tjs) >>= tContract'
  where
@@ -1070,7 +1111,7 @@
 -- utility functions for tensors
 
 nth :: Integer -> [a] -> a
-nth i xs = xs !! (fromIntegral (i - 1))
+nth i xs = xs !! fromIntegral (i - 1)
 
 cdr :: [a] -> [a]
 cdr [] = []
@@ -1080,26 +1121,23 @@
 split w xs
  | V.null xs = []
  | otherwise = let (hs, ts) = V.splitAt (fromIntegral w) xs in
-                 hs:(split w ts)
+                 hs:split w ts
 
 tConcat :: HasTensor a => Index EgisonValue -> [Tensor a] -> EgisonM (Tensor a)
-tConcat s ((Tensor ns@(0:_) _ js):_) = do
-  return $ Tensor (0:ns) V.empty (s:js)
-tConcat s ts@((Tensor ns _ js):_) = return $ Tensor ((fromIntegral (length ts)):ns) (V.concat (map tToVector ts)) (s:js)
+tConcat s (Tensor ns@(0:_) _ js:_) = return $ Tensor (0:ns) V.empty (s:js)
+tConcat s ts@(Tensor ns _ js:_) = return $ Tensor (fromIntegral (length ts):ns) (V.concat (map tToVector ts)) (s:js)
 tConcat s ts = do
   ts' <- mapM getScalar ts
   return $ Tensor [fromIntegral (length ts)] (V.fromList ts') [s]
 
 tConcat' :: HasTensor a => [Tensor a] -> EgisonM (Tensor a)
-tConcat' ((Tensor ns@(0:_) _ _):_) = do
-  return $ Tensor (0:ns) V.empty []
-tConcat' ts@((Tensor ns v _):_) = do
-  return $ Tensor ((fromIntegral (length ts)):ns) (V.concat (map tToVector ts)) []
+tConcat' (Tensor ns@(0:_) _ _:_) = return $ Tensor (0:ns) V.empty []
+tConcat' ts@(Tensor ns v _:_) = return $ Tensor (fromIntegral (length ts):ns) (V.concat (map tToVector ts)) []
 tConcat' ts = do
   ts' <- mapM getScalar ts
   return $ Tensor [fromIntegral (length ts)] (V.fromList ts') []
 
-getScalar :: (Tensor a) -> EgisonM a
+getScalar :: Tensor a -> EgisonM a
 getScalar (Scalar x) = return x
 getScalar _ = throwError $ Default "Inconsitent Tensor order"
 
@@ -1109,7 +1147,7 @@
 findPairs' :: Int -> (a -> a -> Bool) -> [a] -> [(Int, Int)]
 findPairs' _ _ [] = []
 findPairs' m p (x:xs) = case findIndex (p x) xs of
-                    Just i -> (m, m + i + 1):(findPairs' (m + 1) p xs)
+                    Just i -> (m, m + i + 1):findPairs' (m + 1) p xs
                     Nothing -> findPairs' (m + 1) p xs
 
 removePairs :: (Int, Int) -> [a] -> [a]
@@ -1139,14 +1177,16 @@
   show (BoolExpr False) = "#f"
   show (IntegerExpr n) = show n
   show (FloatExpr x y) = showComplexFloat x y
-  show (VarExpr name) = name
+  show (VarExpr name) = show name
   show (PartialVarExpr n) = "%" ++ show n
+  show (FunctionExpr args) = "(function [" ++ unwords (map show args) ++ "])"
+  show (IndexedExpr b expr idxs) = show expr ++ concatMap show idxs
+  show (TupleExpr exprs) = "[" ++ unwords (map show exprs) ++ "]"
 
   show (ApplyExpr fn (TupleExpr [])) = "(" ++ show fn ++ ")"
   show (ApplyExpr fn (TupleExpr args)) = "(" ++ show fn ++ " " ++ unwords (map show args) ++ ")"
   show (ApplyExpr fn arg) = "(" ++ show fn ++ " " ++ show arg ++ ")"
 
-
 instance Show EgisonValue where
   show (Char c) = "c#" ++ [c]
   show (String str) = "\"" ++ T.unpack str ++ "\""
@@ -1154,14 +1194,13 @@
   show (Bool False) = "#f"
   show (ScalarData mExpr) = show mExpr
 --  show (TensorData (Scalar x)) = "invalid scalar:" ++ show x
-  show (TensorData (Tensor [_] xs js)) = "[| " ++ unwords (map show (V.toList xs)) ++ " |]" ++ concat (map show js)
-  show (TensorData (Tensor [0, 0] _ js)) = "[| [|  |] |]" ++ concat (map show js)
-  show (TensorData (Tensor [i, j] xs js)) = "[| " ++ f (fromIntegral j) (V.toList xs) ++ "|]" ++ concat (map show js)
+  show (TensorData (Tensor [_] xs js)) = "[| " ++ unwords (map show (V.toList xs)) ++ " |]" ++ concatMap show js
+  show (TensorData (Tensor [0, 0] _ js)) = "[| [|  |] |]" ++ concatMap show js
+  show (TensorData (Tensor [i, j] xs js)) = "[| " ++ f (fromIntegral j) (V.toList xs) ++ "|]" ++ concatMap show js
    where
     f j [] = ""
     f j xs = "[| " ++ unwords (map show (take j xs)) ++ " |] " ++ f j (drop j xs)
-  show (TensorData (Tensor ns xs js)) = "(tensor {" ++ unwords (map show ns) ++ "} {" ++ unwords (map show (V.toList xs)) ++ "} )" ++ concat (map show js)
-  show (UserIndexedData x js) = show x ++ concat (map show js)
+  show (TensorData (Tensor ns xs js)) = "(tensor {" ++ unwords (map show ns) ++ "} {" ++ unwords (map show (V.toList xs)) ++ "} )" ++ concatMap show js
   show (Float x y) = showComplexFloat x y
   show (InductiveData name []) = "<" ++ name ++ ">"
   show (InductiveData name vals) = "<" ++ name ++ " " ++ unwords (map show vals) ++ ">"
@@ -1176,12 +1215,12 @@
   show (UserMatcher _ BFSMode _) = "#<matcher-bfs>"
   show (UserMatcher _ DFSMode _) = "#<matcher-dfs>"
   show (Func Nothing _ args _) = "(lambda [" ++ unwords (map show args) ++ "] ...)"
-  show (Func (Just name) _ _ _) = name
+  show (Func (Just name) _ _ _) = show name
   show (PartialFunc _ n expr) = show n ++ "#" ++ show expr
   show (CFunc Nothing _ name _) = "(cambda " ++ name ++ " ...)"
-  show (CFunc (Just name) _ _ _) = name
+  show (CFunc (Just name) _ _ _) = show name
   show (MemoizedFunc Nothing _ _ _ names _) = "(memoized-lambda [" ++ unwords names ++ "] ...)"
-  show (MemoizedFunc (Just name) _ _ _ names _) = name
+  show (MemoizedFunc (Just name) _ _ _ names _) = show name
   show (Proc Nothing _ names _) = "(procedure [" ++ unwords names ++ "] ...)"
   show (Proc (Just name) _ _ _) = name
   show (Macro names _) = "(macro [" ++ unwords names ++ "] ...)"
@@ -1201,7 +1240,7 @@
   show (TensorArg name) = "%" ++ name
 
 instance Show ScalarData where
-  show (Div p1 (Plus [(Term 1 [])])) = show p1
+  show (Div p1 (Plus [Term 1 []])) = show p1
   show (Div p1 p2) = "(/ " ++ show p1 ++ " " ++ show p2 ++ ")"
 
 instance Show PolyExpr where
@@ -1222,9 +1261,11 @@
 instance Show SymbolExpr where
   show (Symbol _ (':':':':':':_) []) = "#"
   show (Symbol _ s []) = s
-  show (Symbol _ s js) = s ++ concat (map show js)
+  show (Symbol _ s js) = s ++ concatMap show js
   show (Apply fn mExprs) = "(" ++ show fn ++ " " ++ unwords (map show mExprs) ++ ")"
   show (Quote mExprs) = "'" ++ show mExprs
+  show (FunctionData Nothing argnames args js) = "(function [" ++ unwords (map show argnames) ++ "])" ++ concatMap show js
+  show (FunctionData (Just name) argnames args js) = show name ++ concatMap show js
 
 showComplex :: (Num a, Eq a, Ord a, Show a) => a -> a -> String
 showComplex x 0 = show x
@@ -1247,7 +1288,7 @@
  (Char c) == (Char c') = c == c'
  (String str) == (String str') = str == str'
  (Bool b) == (Bool b') = b == b'
- (ScalarData x) == (ScalarData y) = (x == y)
+ (ScalarData x) == (ScalarData y) = x == y
  (TensorData (Tensor js xs _)) == (TensorData (Tensor js' xs' _)) = (js == js') && (xs == xs')
  (Float x y) == (Float x' y') = (x == x') && (y == y')
  (InductiveData name vals) == (InductiveData name' vals') = (name == name') && (vals == vals')
@@ -1258,7 +1299,6 @@
  (CharHash vals) == (CharHash vals') = vals == vals'
  (StrHash vals) == (StrHash vals') = vals == vals'
  (PrimitiveFunc name1 _) == (PrimitiveFunc name2 _) = name1 == name2
- (UserIndexedData val is) == (UserIndexedData val' is') = (val == val') && (is == is')
  -- Temporary: searching a better solution
  (Func Nothing _ xs1 expr1) == (Func Nothing _ xs2 expr2) = (xs1 == xs2) && (expr1 == expr2)
  (Func (Just name1) _ _ _) == (Func (Just name2) _ _ _) = name1 == name2
@@ -1275,24 +1315,24 @@
   fromEgison :: EgisonValue -> EgisonM a
 
 instance EgisonData Char where
-  toEgison c = Char c
+  toEgison = Char
   fromEgison = liftError . fromCharValue
 
 instance EgisonData Text where
-  toEgison str = String str
+  toEgison = String
   fromEgison = liftError . fromStringValue
 
 instance EgisonData Bool where
-  toEgison b = Bool b
+  toEgison = Bool
   fromEgison = liftError . fromBoolValue
 
 instance EgisonData Integer where
-  toEgison 0 = ScalarData $ mathNormalize' (Div (Plus []) (Plus [(Term 1 [])]))
-  toEgison i = ScalarData $ mathNormalize' (Div (Plus [(Term i [])]) (Plus [(Term 1 [])]))
+  toEgison 0 = ScalarData $ mathNormalize' (Div (Plus []) (Plus [Term 1 []]))
+  toEgison i = ScalarData $ mathNormalize' (Div (Plus [Term i []]) (Plus [Term 1 []]))
   fromEgison = liftError . fromIntegerValue
 
 instance EgisonData Rational where
-  toEgison r = ScalarData $ mathNormalize' (Div (Plus [(Term (numerator r) [])]) (Plus [(Term (denominator r) [])]))
+  toEgison r = ScalarData $ mathNormalize' (Div (Plus [Term (numerator r) []]) (Plus [Term (denominator r) []]))
   fromEgison = liftError . fromRationalValue
 
 instance EgisonData Double where
@@ -1300,7 +1340,7 @@
   fromEgison = liftError . fromFloatValue
 
 instance EgisonData Handle where
-  toEgison h = Port h
+  toEgison = Port
   fromEgison = liftError . fromPortValue
 
 instance (EgisonData a) => EgisonData [a] where
@@ -1315,12 +1355,12 @@
 
 instance (EgisonData a, EgisonData b) => EgisonData (a, b) where
   toEgison (x, y) = Tuple [toEgison x, toEgison y]
-  fromEgison (Tuple (x:y:[])) = (liftM2 (,)) (fromEgison x) (fromEgison y)
+  fromEgison (Tuple [x, y]) = liftM2 (,) (fromEgison x) (fromEgison y)
   fromEgison val = liftError $ throwError $ TypeMismatch "two elements tuple" (Value val)
 
 instance (EgisonData a, EgisonData b, EgisonData c) => EgisonData (a, b, c) where
   toEgison (x, y, z) = Tuple [toEgison x, toEgison y, toEgison z]
-  fromEgison (Tuple (x:y:z:[])) = do
+  fromEgison (Tuple [x, y, z]) = do
     x' <- fromEgison x
     y' <- fromEgison y
     z' <- fromEgison z
@@ -1329,7 +1369,7 @@
 
 instance (EgisonData a, EgisonData b, EgisonData c, EgisonData d) => EgisonData (a, b, c, d) where
   toEgison (x, y, z, w) = Tuple [toEgison x, toEgison y, toEgison z, toEgison w]
-  fromEgison (Tuple (x:y:z:w:[])) = do
+  fromEgison (Tuple [x, y, z, w]) = do
     x' <- fromEgison x
     y' <- fromEgison y
     z' <- fromEgison z
@@ -1350,13 +1390,13 @@
 fromBoolValue val = throwError $ TypeMismatch "bool" (Value val)
 
 fromIntegerValue :: EgisonValue -> Either EgisonError Integer
-fromIntegerValue (ScalarData (Div (Plus []) (Plus [(Term 1 [])]))) = return 0
-fromIntegerValue (ScalarData (Div (Plus [(Term x [])]) (Plus [(Term 1 [])]))) = return x
+fromIntegerValue (ScalarData (Div (Plus []) (Plus [Term 1 []]))) = return 0
+fromIntegerValue (ScalarData (Div (Plus [Term x []]) (Plus [Term 1 []]))) = return x
 fromIntegerValue val = throwError $ TypeMismatch "integer" (Value val)
 
 fromRationalValue :: EgisonValue -> Either EgisonError Rational
 fromRationalValue (ScalarData (Div (Plus []) _)) = return 0
-fromRationalValue (ScalarData (Div (Plus [(Term x [])]) (Plus [(Term y [])]))) = return (x % y)
+fromRationalValue (ScalarData (Div (Plus [Term x []]) (Plus [Term y []]))) = return (x % y)
 fromRationalValue val = throwError $ TypeMismatch "rational" (Value val)
 
 fromFloatValue :: EgisonValue -> Either EgisonError Double
@@ -1454,8 +1494,8 @@
 fromBoolWHNF whnf = throwError $ TypeMismatch "bool" whnf
 
 fromIntegerWHNF :: WHNFData -> Either EgisonError Integer
-fromIntegerWHNF (Value (ScalarData (Div (Plus []) (Plus [(Term 1 [])])))) = return 0
-fromIntegerWHNF (Value (ScalarData (Div (Plus [(Term x [])]) (Plus [(Term 1 [])])))) = return x
+fromIntegerWHNF (Value (ScalarData (Div (Plus []) (Plus [Term 1 []])))) = return 0
+fromIntegerWHNF (Value (ScalarData (Div (Plus [Term x []]) (Plus [Term 1 []])))) = return x
 fromIntegerWHNF whnf = throwError $ TypeMismatch "integer" whnf
 
 fromFloatWHNF :: WHNFData -> Either EgisonError Double
@@ -1474,69 +1514,75 @@
 -- Environment
 --
 
-data Env = Env [HashMap String ObjectRef]
+data Env = Env [HashMap Var ObjectRef] (Maybe VarWithIndices)
  deriving (Show)
 
-data Var = Var String [Index ()]
- deriving (Eq)
-type Binding = (String, ObjectRef)
+data Var = Var [String] [Index ()]
+  deriving (Eq, Generic)
 
-data VarWithIndices = VarWithIndices String [Index String]
+data VarWithIndices = VarWithIndices [String] [Index String]
  deriving (Eq)
 
+instance Hashable (Index ())
+instance Hashable Var
+
+type Binding = (Var, ObjectRef)
+
 instance Show Var where
-  show (Var x is) = x ++ concat (map show is)
+  show (Var xs is) = intercalate "." xs ++ concatMap show is
 
 instance Show VarWithIndices where
-  show (VarWithIndices x is) = x ++ concat (map show is)
+  show (VarWithIndices xs is) = intercalate "." xs ++ concatMap show is
 
 instance Show (Index ()) where
   show (Superscript ()) = "~"
   show (Subscript ()) = "_"
   show (SupSubscript ()) = "~_"
   show (DFscript _ _) = ""
+  show (Userscript _) = "|"
 
 instance Show (Index String) where
   show (Superscript s) = "~" ++ s
   show (Subscript s) = "_" ++ s
   show (SupSubscript s) = "~_" ++ s
   show (DFscript _ _) = ""
+  show (Userscript i) = "|" ++ show i
 
 instance Show (Index EgisonExpr) where
   show (Superscript i) = "~" ++ show i
   show (Subscript i) = "_" ++ show i
   show (SupSubscript i) = "~_" ++ show i
   show (DFscript _ _) = ""
+  show (Userscript i) = "|" ++ show i
 
 instance Show (Index ScalarData) where
   show (Superscript i) = "~" ++ show i
   show (Subscript i) = "_" ++ show i
   show (SupSubscript i) = "~_" ++ show i
   show (DFscript _ _) = ""
+  show (Userscript i) = "|" ++ show i
 
 instance Show (Index EgisonValue) where
-  show (Superscript i) = "~" ++ show i
-  show (Subscript i) = "_" ++ show i
+  show (Superscript i) = case i of
+                         ScalarData (Div (Plus [Term 1 [(Symbol id name (a:indices), 1)]]) (Plus [Term 1 []])) -> "~[" ++ show i ++ "]"
+                         _ -> "~" ++ show i
+  show (Subscript i) = case i of
+                         ScalarData (Div (Plus [Term 1 [(Symbol id name (a:indices), 1)]]) (Plus [Term 1 []])) -> "_[" ++ show i ++ "]"
+                         _ -> "_" ++ show i
   show (SupSubscript i) = "~_" ++ show i
   show (DFscript i j) = "_d" ++ show i ++ show j
-
-instance Show (UserIndex EgisonExpr) where
-  show (Userscript i) = "|" ++ show i
-
-instance Show (UserIndex ScalarData) where
-  show (Userscript i) = "|" ++ show i
-
-instance Show (UserIndex EgisonValue) where
-  show (Userscript i) = "|" ++ show i
+  show (Userscript i) = case i of
+                         ScalarData (Div (Plus [Term 1 [(Symbol id name (a:indices), 1)]]) (Plus [Term 1 []])) -> "_[" ++ show i ++ "]"
+                         _ -> "|" ++ show i
 
 nullEnv :: Env
-nullEnv = Env []
+nullEnv = Env [] Nothing
 
 extendEnv :: Env -> [Binding] -> Env
-extendEnv (Env env) = Env . (: env) . HashMap.fromList
+extendEnv (Env env idx) bdg = Env ((: env) $ HashMap.fromList bdg) idx
 
-refVar :: Env -> String -> Maybe ObjectRef
-refVar (Env env) var = msum $ map (HashMap.lookup var) env
+refVar :: Env -> Var -> Maybe ObjectRef
+refVar (Env env idx) var = msum $ map (HashMap.lookup var) env
 
 --
 -- Pattern Match
@@ -1617,7 +1663,7 @@
 --
 
 newtype EgisonM a = EgisonM {
-    unEgisonM :: (ExceptT EgisonError (FreshT IO) a)
+    unEgisonM :: ExceptT EgisonError (FreshT IO) a
   } deriving (Functor, Applicative, Monad, MonadIO, MonadError EgisonError, MonadFresh)
 
 parallelMapM :: (a -> EgisonM b) -> [a] -> EgisonM [b]
@@ -1629,7 +1675,7 @@
 
 unsafePerformEgison :: (Int, Int) -> EgisonM a -> a
 unsafePerformEgison (x, y) ma =
-  let ((Right ret), _) = unsafePerformIO $ runFreshT (x, y + 1) $ runEgisonM ma in
+  let (Right ret, _) = unsafePerformIO $ runFreshT (x, y + 1) $ runEgisonM ma in
   ret
 --    f' :: (Either EgisonError a) -> (Either EgisonError b) -> EgisonM c
 --    f' (Right x) (Right y) = f x y
@@ -1644,7 +1690,7 @@
   s <- get
   (a, s') <- return $ runFresh s m
   put s'
-  return $ either throwError return $ a   
+  return $ either throwError return a   
   
 fromEgisonM :: EgisonM a -> IO (Either EgisonError a)
 fromEgisonM = modifyCounter . runEgisonM
@@ -1672,11 +1718,13 @@
 
 class (Applicative m, Monad m) => MonadFresh m where
   fresh :: m String
+  freshV :: m Var
 
 instance (Applicative m, Monad m) => MonadFresh (FreshT m) where
   fresh = FreshT $ do (x, y) <- get; modify (\(x,y) -> (x + 1, y))
-                      return $ "$_" ++ (show x) ++ (show y)
-
+                      return $ "$_" ++ show x ++ show y
+  freshV = FreshT $ do (x, y) <- get; modify (\(x,y) -> (x + 1, y))
+                       return $ Var ["$_" ++ show x ++ show y] []
 instance (MonadError e m) => MonadError e (FreshT m) where
   throwError = lift . throwError
   catchError m h = FreshT $ catchError (unFreshT m) (unFreshT . h)
@@ -1687,21 +1735,25 @@
 
 instance (MonadFresh m) => MonadFresh (StateT s m) where
   fresh = lift $ fresh
+  freshV = lift $ freshV
 
 instance (MonadFresh m) => MonadFresh (ExceptT e m) where
   fresh = lift $ fresh
+  freshV = lift $ freshV
 
 instance (MonadFresh m, Monoid e) => MonadFresh (ReaderT e m) where
   fresh = lift $ fresh
+  freshV = lift $ freshV
 
 instance (MonadFresh m, Monoid e) => MonadFresh (WriterT e m) where
   fresh = lift $ fresh
+  freshV = lift $ freshV
 
 instance MonadIO (FreshT IO) where
   liftIO = lift
 
 runFreshT :: Monad m => (Int, Int) -> FreshT m a -> m (a, (Int, Int))
-runFreshT seed = flip (runStateT . unFreshT) seed
+runFreshT = flip (runStateT . unFreshT)
 
 runFresh :: (Int, Int) -> Fresh a -> (a, (Int, Int))
 runFresh seed m = runIdentity $ flip runStateT seed $ unFreshT m
@@ -1728,7 +1780,7 @@
 
 fromMList :: Monad m => MList m a -> m [a]
 fromMList = mfoldr f $ return []
- where f x xs = xs >>= return . (x:)
+  where f x xs = (x:) <$> xs
 
 msingleton :: Monad m => a -> MList m a
 msingleton = flip MCons $ return MNil
@@ -1745,7 +1797,7 @@
 
 mmap :: Monad m => (a -> m b) -> MList m a -> m (MList m b)
 mmap f = mfoldr g $ return MNil
- where g x xs = f x >>= return . flip MCons xs
+  where g x xs = f x >>= return . flip MCons xs
 
 mfor :: Monad m => MList m a -> (a -> m b) -> m (MList m b)
 mfor = flip mmap
@@ -1760,23 +1812,23 @@
 isBool' (Value val) = return $ Value $ Bool $ isBool val
 
 isInteger :: EgisonValue -> Bool
-isInteger (ScalarData (Div (Plus []) (Plus [(Term 1 [])]))) = True
-isInteger (ScalarData (Div (Plus [(Term _ [])]) (Plus [(Term 1 [])]))) = True
+isInteger (ScalarData (Div (Plus []) (Plus [Term 1 []]))) = True
+isInteger (ScalarData (Div (Plus [Term _ []]) (Plus [Term 1 []]))) = True
 isInteger _ = False
 
 isInteger' :: PrimitiveFunc
 isInteger' (Value val) = return $ Value $ Bool $ isInteger val
 
 isRational :: EgisonValue -> Bool
-isRational (ScalarData (Div (Plus []) (Plus [(Term _ [])]))) = True
-isRational (ScalarData (Div (Plus [(Term _ [])]) (Plus [(Term _ [])]))) = True
+isRational (ScalarData (Div (Plus []) (Plus [Term _ []]))) = True
+isRational (ScalarData (Div (Plus [Term _ []]) (Plus [Term _ []]))) = True
 isRational _ = False
 
 isRational' :: PrimitiveFunc
 isRational' (Value val) = return $ Value $ Bool $ isRational val
 
 isSymbol :: EgisonValue -> Bool
-isSymbol (ScalarData (Div (Plus [(Term 1 [(Symbol _ _ _, 1)])]) (Plus [(Term 1 [])]))) = True
+isSymbol (ScalarData (Div (Plus [Term 1 [(Symbol _ _ _, 1)]]) (Plus [Term 1 []]))) = True
 isSymbol _ = False
 
 isScalar :: EgisonValue -> Bool
@@ -1836,8 +1888,19 @@
 isHash' (Intermediate (IStrHash _)) = return $ Value $ Bool True
 isHash' _ = return $ Value $ Bool False
 
-readUTF8File :: FilePath -> (IO String)
+readUTF8File :: FilePath -> IO String
 readUTF8File name = do
   h <- openFile name ReadMode
   hSetEncoding h utf8
   hGetContents h
+
+stringToVar :: String -> Var
+stringToVar name = Var (splitOn "." name) []
+
+varToVarWithIndices :: Var -> VarWithIndices
+varToVarWithIndices (Var xs is) = VarWithIndices xs $ map f is
+ where 
+   f :: Index () -> Index String
+   f (Superscript ()) = Superscript ""
+   f (Subscript ()) = Subscript ""
+   f (SupSubscript ()) = SupSubscript ""
diff --git a/lib/core/base.egi b/lib/core/base.egi
--- a/lib/core/base.egi
+++ b/lib/core/base.egi
@@ -38,7 +38,7 @@
     (lambda $x
       (foldl 2#(%2 %1) x fs))))
 
-(define $flip (lambda [$fn] (lambda [%x %y] (fn y x))))
+(define $flip (lambda [$fn] (lambda [$x $y] (fn y x))))
 
 (define $ref
   (lambda [%xa $is]
diff --git a/lib/core/collection.egi b/lib/core/collection.egi
--- a/lib/core/collection.egi
+++ b/lib/core/collection.egi
@@ -167,7 +167,7 @@
 ;;
 (define $length
   (lambda [$xs]
-    (foldl (lambda [$r $x] (+ r 1)) 0 xs)))
+    (foldl 2#(+ %1 1) 0 xs)))
 
 (define $map
   (lambda [$fn $xs]
@@ -419,14 +419,14 @@
       {@xs x})))
 
 (define $delete-first
-  (lambda [$x $xs]
+  (lambda [%x $xs]
     (match xs (list something)
       {[<nil> {}]
        [<cons ,x $rs> rs]
        [<cons $y $rs> {y @(delete-first x rs)}]})))
 
 (define $delete-first/m
-  (lambda [$a $x $xs]
+  (lambda [$a %x $xs]
     (match xs (list a)
       {[<nil> {}]
        [<cons ,x $rs> rs]
@@ -587,10 +587,15 @@
 ;;
 ;; set operation
 ;;
-(define $unique
+(define $fast-unique
   (lambda [$xs]
     (match-all (sort xs) (list something)
       [<join _ <cons $x !<cons ,x _>>> x])))
+
+(define $unique
+  (lambda [$xs]
+    (reverse (match-all (reverse xs) (list something)
+      [<join _ <cons $x !<join _ <cons ,x _>>>> x]))))
 
 (define $unique/m
   (lambda [$a $xs]
diff --git a/lib/math/algebra/root.egi b/lib/math/algebra/root.egi
--- a/lib/math/algebra/root.egi
+++ b/lib/math/algebra/root.egi
@@ -14,13 +14,13 @@
       (match x math-expr
         {[,0 0]
          [?monomial? (rt-monomial n x)]
-;         [<div <poly $xs> <poly $ys>>
-;          (let {[$xd (reduce gcd xs)]
-;                [$yd (reduce gcd ys)]}
-;            (let {[$d (rt-monomial n (/ xd yd))]}
-;              (*' d
-;                 (rt'' n (*' (/' (sum' (map (/' $ xd) xs)) (sum' (map (/' $ yd) ys)))))
-;                 )))]
+         [<div <poly $xs> <poly $ys>>
+          (let {[$xd (reduce gcd xs)]
+                [$yd (reduce gcd ys)]}
+            (let {[$d (rt-monomial n (/ xd yd))]}
+              (*' d
+                 (rt'' n (*' (/' (sum' (map (/' $ xd) xs)) (sum' (map (/' $ yd) ys)))))
+                 )))]
          [_ (rt'' n x)]})
       (rt'' n x))))
 
diff --git a/lib/math/analysis/derivative.egi b/lib/math/analysis/derivative.egi
--- a/lib/math/analysis/derivative.egi
+++ b/lib/math/analysis/derivative.egi
@@ -10,15 +10,18 @@
       {; symbol
        [,x 1]
        [?symbol? 0]
+       ; function expression
+       [<func _ $argnames $args _> (sum (map2 (lambda [$s $r] (* (user-refs f {s}) (∂/∂ r x))) argnames args))]
        ; function application
        [(,exp $g) (* (exp g) (∂/∂ g x))]
        [(,log $g) (* (/ 1 g) (∂/∂ g x))]
-       [(,cos $g) (* (* -1 (sin g)) (∂/∂ g x))]
-       [(,sin $g) (* (cos g) (∂/∂ g x))]
        [(,sqrt $g) (* (/ 1 (* 2 (sqrt g))) (∂/∂ g x))]
        [(,** $g $h) (* f (∂/∂ (* (log g) h) x))]
+       [(,cos $g) (* (* -1 (sin g)) (∂/∂ g x))]
+       [(,sin $g) (* (cos g) (∂/∂ g x))]
+       [(,arccos $g) (* (/ 1 (sqrt (- 1 (** g 2)))) (∂/∂ g x))]
        [<apply $g $args>
-        (sum (map 2#(* (capply `(add-user-script g %1) args) (∂/∂ %2 x))
+        (sum (map 2#(* (capply `(user-refs g {%1}) args) (∂/∂ %2 x))
                   (zip nats args)))]
        ; quote
        [<quote $g>
diff --git a/lib/math/analysis/integral.egi b/lib/math/analysis/integral.egi
--- a/lib/math/analysis/integral.egi
+++ b/lib/math/analysis/integral.egi
@@ -9,7 +9,7 @@
     (match f math-expr
       {; symbols
        [,x (* (/ 1 2) x^2)]
-       [<symbol _> (* f x)]
+       [<symbol _ _> (* f x)]
        ; function application
        [(,exp ,x) (exp x)]
        [(,cos ,x) (sin x)]
diff --git a/lib/math/expression.egi b/lib/math/expression.egi
--- a/lib/math/expression.egi
+++ b/lib/math/expression.egi
@@ -31,9 +31,9 @@
      [<mult $ $> [integer mult-expr]
       {[<Div <Plus {<Term $n $xs> @{}}> <Plus {<Term 1 {}> @{}}>> {[n (product' (map 2#(**' (to-math-expr' %1) %2) xs))]}]
        [_ {}]}]
-     [<symbol $> [eq]
-      {[<Div <Plus {<Term 1 {[<Symbol $v {}> 1] @{}}> @{}}> <Plus {<Term 1 {}> @{}}>> {v}]
-       [_ {}]}]
+;     [<symbol $> [eq]
+;      {[<Div <Plus {<Term 1 {[<Symbol $v {}> 1] @{}}> @{}}> <Plus {<Term 1 {}> @{}}>> {v}]
+;       [_ {}]}]
      [<symbol $ $> [eq (list index-expr)]
       {[<Div <Plus {<Term 1 {[<Symbol $v $js> 1] @{}}> @{}}> <Plus {<Term 1 {}> @{}}>> {[v js]}]
        [_ {}]}]
@@ -47,13 +47,16 @@
              <Plus {<Term 1 {}> @{}}>>
         {(to-math-expr' mexpr)}]
        [_ {}]}]
+     [<func $ $ $ $> [math-expr (list math-expr) (list math-expr) (list index-expr)]
+      {[<Div <Plus {<Term 1 {[<Function $name $argnames $args $js> 1] @{}}> @{}}> <Plus {<Term 1 {}> @{}}>> {[name argnames args js]}]
+       [_ {}]}]
      [$ [something]
       {[$tgt {(to-math-expr' tgt)}]}]
      }))
 
 (define $index-expr
   (algebraic-data-matcher
-    {<sub math-expr> <sup math-expr>}))
+     {<sub math-expr> <sup math-expr> <user math-expr>}))
 
 (define $poly-expr math-expr)
 (define $term-expr math-expr)
@@ -282,7 +285,7 @@
 (define $find-symbols-from-poly
   (lambda [$poly]
     (match-all poly math-expr
-      [<poly <cons <term _ <cons (& <symbol _> $s) _>> _>> s])))
+      [<poly <cons <term _ <cons (& <symbol _ _> $s) _>> _>> s])))
 
 ;;;
 ;;; Substitute
diff --git a/lib/math/geometry/differential-form.egi b/lib/math/geometry/differential-form.egi
--- a/lib/math/geometry/differential-form.egi
+++ b/lib/math/geometry/differential-form.egi
@@ -14,3 +14,15 @@
 (define $Lie.wedge
   (lambda [%X %Y]
     (- !(. X Y) !(. Y X))))
+
+(define $ι
+  (lambda [%X %Y]
+    (with-symbols {i}
+      (* (df-order Y) (. X...~i (df-normalize Y..._i))))))
+
+(define $Lie
+  (lambda [%X %Y]
+    (match (df-order Y) integer
+      {[,0 (ι X (d Y))]
+       [,N (d (ι X Y))]
+       [_ (+ (ι X (d Y)) (d (ι X Y)))]})))
diff --git a/lib/math/normalize.egi b/lib/math/normalize.egi
--- a/lib/math/normalize.egi
+++ b/lib/math/normalize.egi
@@ -17,14 +17,15 @@
    [rewrite-rule-for-i 1#(contain-symbol? i %1)]
    [rewrite-rule-for-w-term 1#(contain-symbol? w %1)]
    [rewrite-rule-for-rtu-term 1#(contain-function? rtu %1)]
-   [rewrite-rule-for-exp-term 1#(contain-function? exp %1)]
-   [rewrite-rule-for-**-term 1#(contain-function? ** %1)]
+   [rewrite-rule-for-** 1#(contain-function? ** %1)]
+   [rewrite-rule-for-exp 1#(contain-function? exp %1)]
    [rewrite-rule-for-w-poly 1#(contain-symbol? w %1)]
    [rewrite-rule-for-rtu-poly 1#(contain-function? rtu %1)]
    [rewrite-rule-for-sqrt 1#(contain-function? sqrt %1)]
    [rewrite-rule-for-rt 1#(contain-function? rt %1)]
-   ;   [rewrite-rule-for-cos-and-sin 1#(or (contain-function-with-order? cos 2 %1) (contain-function-with-order? sin 2 %1))]
+;   [rewrite-rule-for-cos-and-sin 1#(or (contain-function-with-order? cos 2 %1) (contain-function-with-order? sin 2 %1))]
    [rewrite-rule-for-cos-to-sin 1#(contain-function-with-order? cos 2 %1)]
+   [rewrite-rule-for-d/d 1##t]
    })
 
 ;;
@@ -155,7 +156,7 @@
        [_ term]})))
 
 ;;
-;; exp
+;; **
 ;;
 
 (define $rewrite-rule-for-** (map-terms rewrite-rule-for-**-term $))
@@ -163,7 +164,9 @@
 (define $rewrite-rule-for-**-term
   (lambda [$term]
     (match term math-expr
-      {[(* $a (,** $x $y)^(& ?(gte? $ 2) $n) $r)
+      {[(* $a (,** ,1 _)^_ $r)
+        (rewrite-rule-for-** (*' a r))]
+       [(* $a (,** $x $y)^(& ?(gte? $ 2) $n) $r)
         (rewrite-rule-for-** (*' a (** x (* y n)) r))]
        [(* $a (,** $x $y) (,** ,x $z) $r)
         (rewrite-rule-for-** (*' a (** x (+ y z)) r))]
@@ -253,45 +256,13 @@
 (define $rewrite-rule-for-d/d-poly
   (lambda [$poly]
     (match poly math-expr
-      {[(+ (* $a <apply (& $g <symbol $f $subs>) $args>^$n $mr)
-           (* $b <apply <symbol ,f ?1#(eq?/m (multiset something) subs %1)> ,args>^,n ,mr)
-           $pr)
-        (+ (* (+ a b) (`g args)^n mr)
-           pr)]
-              [_ poly]})))
-
-(define $rewrite-rule-for-cos-to-sin-term'
-  (lambda [$term]
-    (match term math-expr
-      {[(* $a (,cos $x)^,2 $mr)
-        (*' a (-' 1 (sin x)^2) (rewrite-rule-for-cos-to-sin-term' mr))]
-       [_ term]})))
-
-;;
-;; d
-;;
-
-(define $rewrite-rule-for-d (map-terms rewrite-rule-for-d-term $))
-
-(define $rewrite-rule-for-d-term
-  (lambda [$term]
-    (match term math-expr
-      {[(* _ (,d _) (,d _) _)
-        0]
-       [_ term]})))
-
-;;
-;; d/d
-;;
-
-(define $rewrite-rule-for-d/d (map-polys rewrite-rule-for-d/d-poly $))
-
-(define $rewrite-rule-for-d/d-poly
-  (lambda [$poly]
-    (match poly math-expr
-      {[(+ (* $a <apply (& $g <symbol $f $subs>) $args>^$n $mr)
-           (* $b <apply <symbol ,f ?1#(eq?/m (multiset something) subs %1)> ,args>^,n ,mr)
+      {
+       [(+ (* $a (& $f <func $g _ $arg $js>)^$n $mr)
+           (* $b <func ,g _ ,arg ?1#(eq?/m (multiset something) js %1)>^,n ,mr)
            $pr)
-        (+ (* (+ a b) (`g args)^n mr)
-           pr)]
-              [_ poly]})))
+       (rewrite-rule-for-d/d-poly (+' (*' (+ a b) f^n mr) pr))]
+;       [(+ (* $a <apply (& ?scalar? $g <symbol $f $subs>) $args>^$n $mr)
+;           (* $b <apply (& ?scalar? <symbol ,f ?1#(eq?/m (multiset something) subs %1)>) ,args>^,n ,mr)
+;           $pr)
+;       (+ (*' (+ a b) (`g args)^n mr) pr)]
+       [_ poly]})))
diff --git a/sample/math/geometry/covariant-exterior-derivative.egi b/sample/math/geometry/covariant-exterior-derivative.egi
new file mode 100644
--- /dev/null
+++ b/sample/math/geometry/covariant-exterior-derivative.egi
@@ -0,0 +1,59 @@
+;;; Parameters and Metric tensor
+
+(define $x [| θ φ |])
+
+(define $g__ [| [| r^2 0 |] [| 0 (* r^2 (sin θ)^2) |] |])
+(define $g~~ [| [| (/ 1 r^2) 0 |] [| 0 (/ 1 (* r^2 (sin θ)^2)) |] |])
+
+;;; Christoffel symbols
+
+(define $Γ_j_k_l
+  (* (/ 1 2)
+     (+ (∂/∂ g_j_l x_k)
+        (∂/∂ g_j_k x_l)
+        (* -1 (∂/∂ g_k_l x_j)))))
+
+(define $Γ~__ (with-symbols {i} (. g~#~i Γ_i_#_#)))
+
+;;; Riemann curvature tensor
+
+(define $R~i_j_k_l
+  (with-symbols {m}
+    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+       (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
+
+R~#_#_1_1;[| [| 0 0 |] [| 0 0 |] |]~#_#
+R~#_#_1_2;[| [| 0 (sin θ)^2 |] [| -1 0 |] |]~#_#
+R~#_#_2_1;[| [| 0 (* -1 (sin θ)^2) |] [| 1 0 |] |]~#_#
+R~#_#_2_2;[| [| 0 0 |] [| 0 0 |] |]~#_#
+
+;;; Connection form
+
+(define $ω Γ~#_#_#)
+
+;;; Curvature form
+
+(define $wedge
+  (lambda [%X %Y]
+    !(. X Y)))
+
+(define $d
+  (lambda [%A]
+    !((flip ∂/∂) x A)))
+
+(define $D
+  (lambda [%A]
+    (with-symbols {i j}
+      (+ (d A) (wedge ω~i_j A)))))
+
+(define $Ω
+  (with-symbols {i j}
+    (df-normalize (+ (d ω~i_j)
+                     (wedge ω~i_k ω~k_j)))))
+
+Ω~#_#_1_1;[| [| 0 0 |] [| 0 0 |] |]~#_#
+Ω~#_#_1_2;[| [| 0 (/ (sin θ)^2 2) |] [| (/ -1 2) 0 |] |]~#_#
+Ω~#_#_2_1;[| [| 0 (/ (* -1 (sin θ)^2) 2) |] [| (/ 1 2) 0 |] |]~#_#
+Ω~#_#_2_2;[| [| 0 0 |] [| 0 0 |] |]~#_#
+
+
diff --git a/sample/math/geometry/curvature-form.egi b/sample/math/geometry/curvature-form.egi
--- a/sample/math/geometry/curvature-form.egi
+++ b/sample/math/geometry/curvature-form.egi
@@ -9,9 +9,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 (define $Γ~__ (with-symbols {i} (. g~#~i Γ_i_#_#)))
 
@@ -19,7 +19,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_1_1;[| [| 0 0 |] [| 0 0 |] |]~#_#
diff --git a/sample/math/geometry/euler-form-of-S2.egi b/sample/math/geometry/euler-form-of-S2.egi
--- a/sample/math/geometry/euler-form-of-S2.egi
+++ b/sample/math/geometry/euler-form-of-S2.egi
@@ -27,9 +27,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 (define $Γ~__ (with-symbols {i} (. g~#~i Γ_i_#_#)))
 
@@ -56,7 +56,7 @@
     !(. X Y)))
 
 (define $Ω
-  (with-symbols {i j}
+  (with-symbols {i j k}
     (df-normalize (+ (d ω~i_j)
                      (wedge ω~i_k ω~k_j)))))
 Ω~#_#_1_2;[| [| 0 (sin θ) |] [| (* -1 (sin θ)) 0 |] |]~#_#
diff --git a/sample/math/geometry/euler-form-of-T2.egi b/sample/math/geometry/euler-form-of-T2.egi
--- a/sample/math/geometry/euler-form-of-T2.egi
+++ b/sample/math/geometry/euler-form-of-T2.egi
@@ -27,9 +27,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 (define $Γ~__ (with-symbols {i} (. g~#~i Γ_i_#_#)))
 
diff --git a/sample/math/geometry/hodge-E3.egi b/sample/math/geometry/hodge-E3.egi
--- a/sample/math/geometry/hodge-E3.egi
+++ b/sample/math/geometry/hodge-E3.egi
@@ -7,8 +7,8 @@
     (let {[$k (df-order A)]}
       (with-symbols {i j}
         (* (sqrt (abs (M.det g_#_#)))
-           (foldl . (. (subrefs A (map 1#j_%1 (between 1 k)))
-                       (subrefs (ε' N k) (map 1#i_%1 (between 1 N))))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
                   (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
 
 (define $dx [| 1 0 0 |])
diff --git a/sample/math/geometry/hodge-Minkowski.egi b/sample/math/geometry/hodge-Minkowski.egi
--- a/sample/math/geometry/hodge-Minkowski.egi
+++ b/sample/math/geometry/hodge-Minkowski.egi
@@ -7,8 +7,8 @@
     (let {[$k (df-order A)]}
       (with-symbols {i j}
         (* (sqrt (abs (M.det g_#_#)))
-           (foldl . (. (subrefs A (map 1#j_%1 (between 1 k)))
-                       (subrefs (ε' N k) (map 1#i_%1 (between 1 N))))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
                   (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
 
 (define $dt [| 1 0 0 0 |])
diff --git a/sample/math/geometry/hodge-laplacian-E3.egi b/sample/math/geometry/hodge-laplacian-E3.egi
new file mode 100644
--- /dev/null
+++ b/sample/math/geometry/hodge-laplacian-E3.egi
@@ -0,0 +1,53 @@
+;;; Parameters and metrics
+
+(define $N 2)
+
+(define $params [|x y|])
+
+(define $g__ [| [| 1 0 |] [| 0 1 |] |])
+(define $g~~ (M.inverse g_#_#))
+
+;;; Hodge Laplacian
+
+(define $d
+  (lambda [%X]
+    !((flip ∂/∂) params X)))
+
+(define $hodge
+  (lambda [%A]
+    (let {[$k (df-order A)]}
+      (with-symbols {i j}
+        (* (sqrt (abs (M.det g_#_#)))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
+                  (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
+
+(define $δ
+  (lambda [%A]
+    (let {[$k (df-order A)]}
+      (* (** -1 (+ (* N (+ k 1)) 1))
+         (hodge (d (hodge A)))))))
+
+(define $Δ
+  (lambda [%A]
+    (match (df-order A) integer
+      {[,0 (δ (d A))]
+       [,2 (d (δ A))]
+       [_ (+ (d (δ A)) (δ (d A)))]})))
+
+(define $f (function [x y]))
+
+(d f)
+;[| f|x f|y |]
+
+(hodge (d f))
+;[| (* -1 f|y) f|x |]
+
+(d (hodge (d f)))
+;[| [| (* -1 f|y|x) f|x|x |] [| (* -1 f|y|y) f|x|y |] |]
+
+(hodge (d (hodge (d f))))
+;(+ f|y|y f|x|x)
+
+(Δ f)
+;(+ (* -1 f|y|y) (* -1 f|x|x))
diff --git a/sample/math/geometry/hodge-laplacian-one-form.egi b/sample/math/geometry/hodge-laplacian-one-form.egi
new file mode 100644
--- /dev/null
+++ b/sample/math/geometry/hodge-laplacian-one-form.egi
@@ -0,0 +1,52 @@
+;;; Parameters and metrics
+
+(define $N 3)
+
+(define $params [| x y z |])
+
+(define $g__ [| [| 1 0 0 |] [| 0 1 0 |] [| 0 0 1 |] |])
+(define $g~~ (M.inverse g_#_#))
+
+;;; Hodge Laplacian
+
+(define $d
+  (lambda [%X]
+    !((flip ∂/∂) params X)))
+
+(define $hodge
+  (lambda [%A]
+    (let {[$k (df-order A)]}
+      (with-symbols {i j}
+        (* (sqrt (abs (M.det g_#_#)))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
+                  (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
+
+(define $δ
+  (lambda [%A]
+    (let {[$r (df-order A)]}
+      (* (** -1 (+ (* N r) 1))
+         (hodge (d (hodge A)))))))
+
+(define $Δ
+  (lambda [%A]
+    (match (df-order A) integer
+      {[,0 (δ (d A))]
+       [,3 (d (δ A))]
+       [_ (+ (d (δ A)) (δ (d A)))]})))
+
+(define $ux (function [t x y z]))
+(define $uy (function [t x y z]))
+(define $uz (function [t x y z]))
+(define $u [| ux uy uz |])
+
+(Δ u)
+;[| (+ ux|x|x ux|z|z ux|y|y) (+ uy|y|y uy|z|z uy|x|x) (+ uz|z|z uz|y|y uz|x|x) |]
+
+(define $vx (function [t x y z]))
+(define $vy (function [t x y z]))
+(define $vz (function [t x y z]))
+(define $v [|[| 0 vz (* -1 vy) |] [| (* -1 vz) 0 vx |] [| vy (* -1 vx) 0 |]|])
+
+(df-normalize (Δ v))
+;[| [| 0 (+ vz|x|x vz|z|z vz|y|y) (+ (* -1 vy|x|x) (* -1 vy|y|y) (* -1 vy|z|z)) |] [| (+ (* -1 vz|y|y) (* -1 vz|x|x) (* -1 vz|z|z)) 0 (+ vx|y|y vx|x|x vx|z|z) |] [| (+ vy|z|z vy|x|x vy|y|y) (+ (* -1 vx|z|z) (* -1 vx|y|y) (* -1 vx|x|x)) 0 |] |]
diff --git a/sample/math/geometry/hodge-laplacian-polar.egi b/sample/math/geometry/hodge-laplacian-polar.egi
--- a/sample/math/geometry/hodge-laplacian-polar.egi
+++ b/sample/math/geometry/hodge-laplacian-polar.egi
@@ -18,14 +18,14 @@
     (let {[$k (df-order A)]}
       (with-symbols {i j}
         (* (sqrt (abs (M.det g_#_#)))
-           (foldl . (. (subrefs A (map 1#j_%1 (between 1 k)))
-                       (subrefs (ε' N k) (map 1#i_%1 (between 1 N))))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
                   (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
 
 (define $δ
   (lambda [%A]
-    (let {[$r (df-order A)]}
-      (* (** -1 (+ (* N r) 1))
+    (let {[$k (df-order A)]}
+      (* (** -1 (+ (* N (+ k 1)) 1))
          (hodge (d (hodge A)))))))
 
 (define $Δ
@@ -35,5 +35,19 @@
        [,2 (d (δ A))]
        [_ (+ (d (δ A)) (δ (d A)))]})))
 
-(Δ (f r θ))
-;(/ (+ (* -1 (f|2|2 r θ)) (* -1 r (f|1 r θ)) (* -1 r^2 (f|1|1 r θ))) r^2)
+(define $f (function [r θ]))
+
+(d f)
+;[| f|r f|θ |]
+
+(hodge (d f))
+;[| (/ (* -1 f|θ) r) (* r f|r) |]
+
+(d (hodge (d f)))
+;[| [| (/ (+ (* -1 f|θ|r r) f|θ) r^2) (+ f|r (* r f|r|r)) |] [| (/ (* -1 f|θ|θ) r) (* r f|r|θ) |] |]
+
+(hodge (d (hodge (d f))))
+;(/ (+ f|θ|θ (* r f|r) (* r^2 f|r|r)) r^2)
+
+(Δ f)
+;(/ (+ (* -1 f|θ|θ) (* -1 r f|r) (* -1 r^2 f|r|r)) r^2)
diff --git a/sample/math/geometry/hodge-laplacian-spherical.egi b/sample/math/geometry/hodge-laplacian-spherical.egi
--- a/sample/math/geometry/hodge-laplacian-spherical.egi
+++ b/sample/math/geometry/hodge-laplacian-spherical.egi
@@ -18,8 +18,8 @@
     (let {[$k (df-order A)]}
       (with-symbols {i j}
         (* (sqrt (abs (M.det g_#_#)))
-           (foldl . (. (subrefs A (map 1#j_%1 (between 1 k)))
-                       (subrefs (ε' N k) (map 1#i_%1 (between 1 N))))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
                   (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
 
 (define $δ
diff --git a/sample/math/geometry/hodge-laplacian.egi b/sample/math/geometry/hodge-laplacian.egi
--- a/sample/math/geometry/hodge-laplacian.egi
+++ b/sample/math/geometry/hodge-laplacian.egi
@@ -18,8 +18,8 @@
     (let {[$k (df-order A)]}
       (with-symbols {i j}
         (* (sqrt (abs (M.det g_#_#)))
-           (foldl . (. (subrefs A (map 1#j_%1 (between 1 k)))
-                       (subrefs (ε' N k) (map 1#i_%1 (between 1 N))))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
                   (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
 
 (define $δ
diff --git a/sample/math/geometry/lie.egi b/sample/math/geometry/lie.egi
new file mode 100644
--- /dev/null
+++ b/sample/math/geometry/lie.egi
@@ -0,0 +1,48 @@
+(define $N 3)
+(define $params [| x y z |])
+(define $g [| [| 1 0 0 |] [| 0 1 0 |] [| 0 0 1 |] |])
+
+(define $d
+  (lambda [%X]
+    !((flip ∂/∂) params X)))
+
+(define $hodge
+  (lambda [%A]
+    (let {[$k (df-order A)]}
+      (with-symbols {i j}
+        (* (sqrt (abs (M.det g_#_#)))
+           (foldl . (. A_[j_1]..._[j_k]
+                       (ε' N k)_[i_1]..._[i_N])
+                  (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
+
+(define $dx [| 1 0 0 |])
+(define $dy [| 0 1 0 |])
+(define $dz [| 0 0 1 |])
+
+(define $ι
+  (lambda [%X %Y]
+    (with-symbols {i}
+      (* (df-order Y) (. X...~i (df-normalize Y..._i))))))
+
+(define $Lie
+  (lambda [%X %Y]
+    (match (df-order Y) integer
+      {[,0 (ι X (d Y))]
+       [,N (d (ι X Y))]
+       [_ (+ (ι X (d Y)) (d (ι X Y)))]})))
+
+(define $ρ (function [t x y z]))
+(define $*ρ (df-normalize (hodge ρ)))
+
+(define $u_ (generate-tensor 1#(function [t x y z]) {3}))
+(define $u [| u_1 u_2 u_3 |])
+
+(df-normalize (+ (∂/∂ *ρ t) (Lie u *ρ)))
+;(tensor {3 3 3} {0 0 0 0 0 (/ (+ ρ|t (* u_1|x ρ) (* u_1 ρ|x) (* u_2|y ρ) (* u_2 ρ|y) (* u_3|z ρ) (* u_3 ρ|z)) 6) 0 (/ (+ (* -1 ρ|t) (* -1 u_1|x ρ) (* -1 u_1 ρ|x) (* -1 u_3|z ρ) (* -1 u_3 ρ|z) (* -1 u_2|y ρ) (* -1 u_2 ρ|y)) 6) 0 0 0 (/ (+ (* -1 ρ|t) (* -1 u_2|y ρ) (* -1 u_2 ρ|y) (* -1 u_1|x ρ) (* -1 u_1 ρ|x) (* -1 u_3|z ρ) (* -1 u_3 ρ|z)) 6) 0 0 0 (/ (+ ρ|t (* u_2|y ρ) (* u_2 ρ|y) (* u_3|z ρ) (* u_3 ρ|z) (* u_1|x ρ) (* u_1 ρ|x)) 6) 0 0 0 (/ (+ ρ|t (* u_3|z ρ) (* u_3 ρ|z) (* u_1|x ρ) (* u_1 ρ|x) (* u_2|y ρ) (* u_2 ρ|y)) 6) 0 (/ (+ (* -1 ρ|t) (* -1 u_3|z ρ) (* -1 u_3 ρ|z) (* -1 u_2|y ρ) (* -1 u_2 ρ|y) (* -1 u_1|x ρ) (* -1 u_1 ρ|x)) 6) 0 0 0 0 0} )
+
+(df-normalize (+ (∂/∂ *ρ t) (Lie u *ρ)))_1_2_3
+;(/ (+ ρ|t
+;      (* u_1|x ρ) (* u_1 ρ|x)
+;      (* u_2|y ρ) (* u_2 ρ|y)
+;      (* u_3|z ρ) (* u_3 ρ|z))
+;   6)
diff --git a/sample/math/geometry/polar-laplacian-2d-2.egi b/sample/math/geometry/polar-laplacian-2d-2.egi
--- a/sample/math/geometry/polar-laplacian-2d-2.egi
+++ b/sample/math/geometry/polar-laplacian-2d-2.egi
@@ -12,7 +12,7 @@
 ;; Local coordinates
 ;;
 
-(define $e ((∂/∂ X~# $) x_#))
+(define $e ((∂/∂ X_# $) x~#))
 e
 ;[| [| (cos θ) (sin θ) |] [| (* -1 r (sin θ)) (* r (cos θ)) |] |]
 
@@ -33,9 +33,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_l x_k)
-          (∂/∂ g_j_k x_l)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_l x~k)
+          (∂/∂ g_j_k x~l)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_#_#_#;(tensor {2 2 2} {0 0 0 (* -1 r) 0 r r 0} )_#_#_#
 Γ_1_#_#;[| [| 0 0 |] [| 0 (* -1 r) |] |]_#_#
@@ -57,12 +57,12 @@
 ;; Derive Laplacian
 ;;
 
-(. g~i~j (∂/∂ (∂/∂ (f r θ) x_j) x_i))
+(. g~i~j (∂/∂ (∂/∂ (f r θ) x~j) x~i))
 ;(/ (+ (* (f|1|1 r θ) r^2) (f|2|2 r θ)) r^2)
-(. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ) x_k))
+(. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ) x~k))
 ;(/ (* -1 (f|1 r θ)) r)
 
-(define $Laplacian (- (. g~i~j (∂/∂ (∂/∂ (f r θ) x_j) x_i))
-                        (. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ) x_k))))
+(define $Laplacian (- (. g~i~j (∂/∂ (∂/∂ (f r θ) x~j) x~i))
+                        (. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ) x~k))))
 Laplacian
 ;(/ (+ (* (f|1|1 r θ) r^2) (f|2|2 r θ) (* (f|1 r θ) r)) r^2)
diff --git a/sample/math/geometry/polar-laplacian-2d-3.egi b/sample/math/geometry/polar-laplacian-2d-3.egi
--- a/sample/math/geometry/polar-laplacian-2d-3.egi
+++ b/sample/math/geometry/polar-laplacian-2d-3.egi
@@ -12,7 +12,7 @@
 ;; Local coordinates
 ;;
 
-(define $e ((∂/∂ X~# $) x_#))
+(define $e ((∂/∂ X_# $) x~#))
 e
 ;[| [| (cos θ) (sin θ) |] [| (* -1 r (sin θ)) (* r (cos θ)) |] |]
 
@@ -33,6 +33,6 @@
 (define $sqrt-g (sqrt (M.det g_#_#)))
 sqrt-g;r
 
-(define $Laplacian (/ (contract + (∂/∂ (* sqrt-g (. g~i~j (∂/∂ (f r θ) x_j))) x_i)) sqrt-g))
+(define $Laplacian (/ (contract + (∂/∂ (* sqrt-g (. g~i~j (∂/∂ (f r θ) x~j))) x~i)) sqrt-g))
 Laplacian
 ;(/ (+ (* (f|1 r θ) r) (* r^2 (f|1|1 r θ)) (f|2|2 r θ)) r^2)
diff --git a/sample/math/geometry/polar-laplacian-3d-2.egi b/sample/math/geometry/polar-laplacian-3d-2.egi
--- a/sample/math/geometry/polar-laplacian-3d-2.egi
+++ b/sample/math/geometry/polar-laplacian-3d-2.egi
@@ -13,7 +13,7 @@
 ;; Local coordinates
 ;;
 
-(define $e ((∂/∂ X~# $) x_#))
+(define $e ((∂/∂ X_# $) x~#))
 e
 ;[|[| (* (sin θ) (cos φ)) (* (sin θ) (sin φ)) (cos θ) |]
 ;  [| (* r (cos θ) (cos φ)) (* r (cos θ) (sin φ)) (* -1 r (sin θ)) |]
@@ -36,9 +36,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_l x_k)
-          (∂/∂ g_j_k x_l)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_l x~k)
+          (∂/∂ g_j_k x~l)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_#_#_#;(tensor {3 3 3} {0 0 0 0 (* -1 r) 0 0 0 (* -1 r (sin θ)^2) 0 r 0 r 0 0 0 0 (* -1 r^2 (sin θ) (cos θ)) 0 0 (* r (sin θ)^2) 0 0 (* r^2 (sin θ) (cos θ)) (* r (sin θ)^2) (* r^2 (sin θ) (cos θ)) 0} )_#_#_#
 Γ_1_#_#;[| [| 0 0 0 |] [| 0 (* -1 r) 0 |] [| 0 0 (* -1 r (sin θ)^2) |] |]_#_#
@@ -62,12 +62,12 @@
 ;; Laplacian
 ;;
 
-(. g~i~j (∂/∂ (∂/∂ (f r θ φ) x_j) x_i))
+(. g~i~j (∂/∂ (∂/∂ (f r θ φ) x~j) x~i))
 ;(/ (+ (* (f|1|1 r θ φ) r^2 (sin θ)^2) (* (f|2|2 r θ φ) (sin θ)^2) (f|3|3 r θ φ)) (* r^2 (sin θ)^2))
-(. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ φ) x_k))
+(. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ φ) x~k))
 ;(/ (+ (* -2 (f|1 r θ φ) r (sin θ)) (* -1 (cos θ) (f|2 r θ φ))) (* r^2 (sin θ)))
 
-(define $Laplacian (- (. g~i~j (∂/∂ (∂/∂ (f r θ φ) x_j) x_i))
-                        (. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ φ) x_k))))
+(define $Laplacian (- (. g~i~j (∂/∂ (∂/∂ (f r θ φ) x~j) x~i))
+                        (. (. g~i~j Γ~k_i_j) (∂/∂ (f r θ φ) x~k))))
 Laplacian
 ;(/ (+ (* (f|1|1 r θ φ) r^2 (sin θ)^2) (* (f|2|2 r θ φ) (sin θ)^2) (f|3|3 r θ φ) (* 2 (f|1 r θ φ) r (sin θ)^2) (* (cos θ) (f|2 r θ φ) (sin θ))) (* r^2 (sin θ)^2))
diff --git a/sample/math/geometry/polar-laplacian-3d-3.egi b/sample/math/geometry/polar-laplacian-3d-3.egi
--- a/sample/math/geometry/polar-laplacian-3d-3.egi
+++ b/sample/math/geometry/polar-laplacian-3d-3.egi
@@ -13,7 +13,7 @@
 ;; Local coordinates
 ;;
 
-(define $e ((∂/∂ X~# $) x_#))
+(define $e ((∂/∂ X_# $) x~#))
 e
 ;[|[| (* (sin θ) (cos φ)) (* (sin θ) (sin φ)) (cos θ) |]
 ;  [| (* r (cos θ) (cos φ)) (* r (cos θ) (sin φ)) (* -1 r (sin θ)) |]
@@ -36,6 +36,6 @@
 (define $sqrt-g (sqrt (M.det g_#_#)))
 sqrt-g;(* r^2 (sin θ))
 
-(define $Laplacian (/ (contract + (∂/∂ (* sqrt-g (. g~i~j (∂/∂ (f r θ φ) x_j))) x_i)) sqrt-g))
+(define $Laplacian (/ (contract + (∂/∂ (* sqrt-g (. g~i~j (∂/∂ (f r θ φ) x~j))) x~i)) sqrt-g))
 Laplacian
 ;(/ (+ (* 2 r (sin θ)^2 (f|1 r θ φ)) (* r^2 (sin θ)^2 (f|1|1 r θ φ)) (* (cos θ) (f|2 r θ φ) (sin θ)) (* (sin θ)^2 (f|2|2 r θ φ)) (f|3|3 r θ φ)) (* (sin θ)^2 r^2))
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-FLRW-metric.egi b/sample/math/geometry/riemann-curvature-tensor-of-FLRW-metric.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-FLRW-metric.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-FLRW-metric.egi
@@ -33,9 +33,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_k x_l)
-        (∂/∂ g_j_l x_k)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_k x~l)
+        (∂/∂ g_j_l x~k)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_1_#_#;[| [| 0 0 0 0 |] [| 0 (/ (* -1 (a w) (a|1 w)) '(+ 1 (* -1 K r^2))) 0 0 |] [| 0 0 (* -1 (a w) (a|1 w) r^2) 0 |] [| 0 0 0 (* -1 (a w) (a|1 w) r^2 (sin θ)^2) |] |]_#_#
 Γ_2_#_#;[| [| 0 (/ (* (a w) (a|1 w)) '(+ 1 (* -1 K r^2))) 0 0 |] [| (/ (* (a w) (a|1 w)) '(+ 1 (* -1 K r^2))) (/ (* K r (a w)^2) '(+ 1 (* -1 K r^2))^2) 0 0 |] [| 0 0 (* -1 (a w)^2 r) 0 |] [| 0 0 0 (* -1 (a w)^2 r (sin θ)^2) |] |]_#_#
@@ -59,7 +59,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_1_1;[| [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] |]~#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-M3-conformal.egi b/sample/math/geometry/riemann-curvature-tensor-of-M3-conformal.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-M3-conformal.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-M3-conformal.egi
@@ -20,9 +20,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_l x_k)
-          (∂/∂ g_j_k x_l)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_l x~k)
+          (∂/∂ g_j_k x~l)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_#_#_#
 
@@ -42,7 +42,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-M5-conformal.egi b/sample/math/geometry/riemann-curvature-tensor-of-M5-conformal.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-M5-conformal.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-M5-conformal.egi
@@ -20,9 +20,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_l x_k)
-          (∂/∂ g_j_k x_l)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_l x~k)
+          (∂/∂ g_j_k x~l)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_#_#_#
 
@@ -42,7 +42,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S1.egi b/sample/math/geometry/riemann-curvature-tensor-of-S1.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S1.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S1.egi
@@ -32,9 +32,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_k x_l)
-          (∂/∂ g_j_l x_k)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_k x~l)
+          (∂/∂ g_j_l x~k)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_#_#_#;(tensor {1 1 1} {0} )_#_#_#
 
@@ -54,7 +54,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#;(tensor {1 1 1 1} {0} )~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S2.egi b/sample/math/geometry/riemann-curvature-tensor-of-S2.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S2.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S2.egi
@@ -13,7 +13,7 @@
 ;; Local basis
 ;;
 
-(define $e_i_j (∂/∂ X_j x_i))
+(define $e_i_j (∂/∂ X_j x~i))
 e_i_j
 ;[|[|(* r (cos θ) (cos φ)) (* r (cos θ) (sin φ)) (* -1 r (sin θ)) |]
 ;  [|(* -1 r (sin θ) (sin φ)) (* r (sin θ) (cos φ)) 0 |]
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-conformal-fast.egi b/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-conformal-fast.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-conformal-fast.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-conformal-fast.egi
@@ -26,9 +26,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_#_#_#
 
@@ -46,7 +46,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-fast.egi b/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-fast.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-fast.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S2xS3-fast.egi
@@ -25,9 +25,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_#_#_#
 
@@ -45,7 +45,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S2xS3.egi b/sample/math/geometry/riemann-curvature-tensor-of-S2xS3.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S2xS3.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S2xS3.egi
@@ -25,9 +25,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_#_#_#
 
@@ -45,7 +45,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S3.egi b/sample/math/geometry/riemann-curvature-tensor-of-S3.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S3.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S3.egi
@@ -37,9 +37,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_k x_l)
-        (∂/∂ g_j_l x_k)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_k x~l)
+        (∂/∂ g_j_l x~k)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_1_#_#;[| [| 0 0 0 |] [| 0 (* -1 r^2 (sin θ) (cos θ)) 0 |] [| 0 0 (* -1 r^2 (sin θ) (cos θ) (sin φ)^2) |] |]_#_#
 Γ_2_#_#;[| [| 0 (* r^2 (sin θ) (cos θ)) 0 |] [| (* r^2 (sin θ) (cos θ)) 0 0 |] [| 0 0 (* -1 r^2 (sin θ)^2 (sin φ) (cos φ)) |] |]_#_#
@@ -61,7 +61,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_1_1;[| [| 0 0 0 |] [| 0 0 0 |] [| 0 0 0 |] |]~#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S4.egi b/sample/math/geometry/riemann-curvature-tensor-of-S4.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S4.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S4.egi
@@ -40,9 +40,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_k x_l)
-        (∂/∂ g_j_l x_k)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_k x~l)
+        (∂/∂ g_j_l x~k)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_1_#_#;[| [| 0 0 0 0 |] [| 0 (/ (* -1 r^2 (sin (* 2 θ))) 2) 0 0 |] [| 0 0 (/ (* -1 r^2 (sin (* 2 θ)) (sin φ)^2) 2) 0 |] [| 0 0 0 (/ (* -1 r^2 (sin (* 2 θ)) (sin φ)^2 (sin ψ)^2) 2) |] |]_#_#
 Γ_2_#_#;[| [| 0 (/ (* r^2 (sin (* 2 θ))) 2) 0 0 |] [| (/ (* r^2 (sin (* 2 θ))) 2) 0 0 0 |] [| 0 0 (/ (* -1 r^2 (sin θ)^2 (sin (* 2 φ))) 2) 0 |] [| 0 0 0 (/ (* -1 r^2 (sin θ)^2 (sin (* 2 φ)) (sin ψ)^2) 2) |] |]_#_#
@@ -66,7 +66,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_1_1;[| [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] |]~#_#
@@ -112,7 +112,7 @@
 
 (define $∇Ric___
   (with-symbols {i j k l m n}
-    (- (∂/∂ Ric_i_j x_m)
+    (- (∂/∂ Ric_i_j x~m)
        (. Γ~n_m_i Ric_n_j)
        (. Γ~n_m_j Ric_i_n)
        )))
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal-weyl.egi b/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal-weyl.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal-weyl.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal-weyl.egi
@@ -42,9 +42,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_#_#_#
 
@@ -62,7 +62,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal.egi b/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S5-conformal.egi
@@ -42,9 +42,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_#_#_#
 
@@ -62,7 +62,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S5-weyl.egi b/sample/math/geometry/riemann-curvature-tensor-of-S5-weyl.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S5-weyl.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S5-weyl.egi
@@ -43,9 +43,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_l x_k)
-          (∂/∂ g_j_k x_l)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_l x~k)
+          (∂/∂ g_j_k x~l)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_#_#_#
 
@@ -65,7 +65,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S5.egi b/sample/math/geometry/riemann-curvature-tensor-of-S5.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S5.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S5.egi
@@ -42,9 +42,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_#_#_#
 
@@ -62,7 +62,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S7-conformal.egi b/sample/math/geometry/riemann-curvature-tensor-of-S7-conformal.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S7-conformal.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S7-conformal.egi
@@ -36,9 +36,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 ;;
 ;; Christoffel symbols of the second kind
@@ -52,7 +52,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 ;;
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-S7.egi b/sample/math/geometry/riemann-curvature-tensor-of-S7.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-S7.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-S7.egi
@@ -39,9 +39,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_l x_k)
-        (∂/∂ g_j_k x_l)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_l x~k)
+        (∂/∂ g_j_k x~l)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 ;;
 ;; Christoffel symbols of the second kind
@@ -55,7 +55,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 ;;
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-Schwarzschild-metric.egi b/sample/math/geometry/riemann-curvature-tensor-of-Schwarzschild-metric.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-Schwarzschild-metric.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-Schwarzschild-metric.egi
@@ -32,9 +32,9 @@
 
 (define $Γ_j_k_l
   (* (/ 1 2)
-     (+ (∂/∂ g_j_k x_l)
-        (∂/∂ g_j_l x_k)
-        (* -1 (∂/∂ g_k_l x_j)))))
+     (+ (∂/∂ g_j_k x~l)
+        (∂/∂ g_j_l x~k)
+        (* -1 (∂/∂ g_k_l x~j)))))
 
 Γ_1_#_#;[| [| 0 (/ (+ (* c^2 r) (* -1 '(+ (* c^2 r) (* -2 G M)))) (* 2 c^2 r^2)) 0 0 |] [| (/ (+ (* c^2 r) (* -1 '(+ (* c^2 r) (* -2 G M)))) (* 2 c^2 r^2)) 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] |]_#_#
 Γ_2_#_#;[| [| (/ (+ (* -1 c^2 r) '(+ (* c^2 r) (* -2 G M))) (* 2 c^2 r^2)) 0 0 0 |] [| 0 (/ (+ (* -1 c^2 '(+ (* c^2 r) (* -2 G M))) (* c^4 r)) (* 2 '(+ (* c^2 r) (* -2 G M))^2)) 0 0 |] [| 0 0 r 0 |] [| 0 0 0 (* r (sin θ)^2) |] |]_#_#
@@ -58,7 +58,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (expand-all (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (expand-all (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
                    (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l))))))
 
 R~#_#_1_1;[| [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] |]~#_#
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-T2.egi b/sample/math/geometry/riemann-curvature-tensor-of-T2.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-T2.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-T2.egi
@@ -35,9 +35,9 @@
 
 (define $Γ_i_j_k
   (* (/ 1 2)
-     (+ (∂/∂ g_i_j x_k)
-        (∂/∂ g_i_k x_j)
-        (* -1 (∂/∂ g_j_k x_i)))))
+     (+ (∂/∂ g_i_j x~k)
+        (∂/∂ g_i_k x~j)
+        (* -1 (∂/∂ g_j_k x~i)))))
 
 Γ_#_#_#;(tensor {2 2 2} {0 0 0 (* '(+ (* a (cos θ)) b) a (sin θ)) 0 (* -1 '(+ (* a (cos θ)) b) a (sin θ)) (* -1 '(+ (* a (cos θ)) b) a (sin θ)) 0} )_#_#_#
 Γ_1_#_#;[| [| 0 0 |] [| 0 (* '(+ (* a (cos θ)) b) a (sin θ)) |] |]_#_#
@@ -58,7 +58,7 @@
 ;;
 (define $∇g___
   (with-symbols {i j m n}
-    (- (∂/∂ g_i_j x_m)
+    (- (∂/∂ g_i_j x~m)
        (. Γ~n_m_i g_n_j)
        (. Γ~n_m_j g_i_n))))
 
@@ -70,7 +70,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#;(tensor {2 2 2 2} {0 0 0 0 0 (/ (* '(+ (* a (cos θ)) b) (cos θ)) a) (/ (* -1 '(+ (* a (cos θ)) b) (cos θ)) a) 0 0 (/ (* -1 a (cos θ)) '(+ (* a (cos θ)) b)) (/ (* a (cos θ)) '(+ (* a (cos θ)) b)) 0 0 0 0 0} )~#_#_#_#
@@ -109,7 +109,7 @@
 
 (define $∇R_____
   (with-symbols {i j k l m n}
-    (- (∂/∂ R_i_j_k_l x_m)
+    (- (∂/∂ R_i_j_k_l x~m)
        (. Γ~n_m_i R_n_j_k_l)
        (. Γ~n_m_j R_i_n_k_l)
        (. Γ~n_m_k R_i_j_n_l)
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-empty-Schwarzschild-spacetime.egi b/sample/math/geometry/riemann-curvature-tensor-of-empty-Schwarzschild-spacetime.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-empty-Schwarzschild-spacetime.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-empty-Schwarzschild-spacetime.egi
@@ -33,9 +33,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_l x_k)
-          (∂/∂ g_j_k x_l)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_l x~k)
+          (∂/∂ g_j_k x~l)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_1_#_#;[| [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 0 0 |] |]_#_#
 Γ_2_#_#;[| [| 0 0 0 0 |] [| 0 0 0 0 |] [| 0 0 r 0 |] [| 0 0 0 (* r (sin θ)^2) |] |]_#_#
@@ -61,7 +61,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 
diff --git a/sample/math/geometry/riemann-curvature-tensor-of-spherical-space.egi b/sample/math/geometry/riemann-curvature-tensor-of-spherical-space.egi
--- a/sample/math/geometry/riemann-curvature-tensor-of-spherical-space.egi
+++ b/sample/math/geometry/riemann-curvature-tensor-of-spherical-space.egi
@@ -30,9 +30,9 @@
 (define $Γ___
   (with-symbols {j k l}
     (* (/ 1 2)
-       (+ (∂/∂ g_j_l x_k)
-          (∂/∂ g_j_k x_l)
-          (* -1 (∂/∂ g_k_l x_j))))))
+       (+ (∂/∂ g_j_l x~k)
+          (∂/∂ g_j_k x~l)
+          (* -1 (∂/∂ g_k_l x~j))))))
 
 Γ_1_#_#;[| [| 0 0 0 |] [| 0 (* -1 r) 0 |] [| 0 0 (* -1 r (sin θ)^2) |] |]_#_#
 Γ_2_#_#;[| [| 0 r 0 |] [| r 0 0 |] [| 0 0 (* -1 r^2 (sin θ) (cos θ)) |] |]_#_#
@@ -56,7 +56,7 @@
 
 (define $R~i_j_k_l
   (with-symbols {m}
-    (+ (- (∂/∂ Γ~i_j_l x_k) (∂/∂ Γ~i_j_k x_l))
+    (+ (- (∂/∂ Γ~i_j_l x~k) (∂/∂ Γ~i_j_k x~l))
        (- (. Γ~m_j_l Γ~i_m_k) (. Γ~m_j_k Γ~i_m_l)))))
 
 R~#_#_#_#;(tensor {3 3 3 3} {0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0} )~#_#_#_#
diff --git a/sample/math/geometry/yang-mills-equation-of-U1-gauge-theory.egi b/sample/math/geometry/yang-mills-equation-of-U1-gauge-theory.egi
--- a/sample/math/geometry/yang-mills-equation-of-U1-gauge-theory.egi
+++ b/sample/math/geometry/yang-mills-equation-of-U1-gauge-theory.egi
@@ -11,8 +11,8 @@
     (let {[$k (df-order A)]}
       (with-symbols {i j}
         (* (sqrt (abs (M.det g_#_#)))
-           (foldl . (. (subrefs A (map 1#j_%1 (between 1 k)))
-                       (subrefs (ε' N k) (map 1#i_%1 (between 1 N))))
+           (foldl . (. (ε' N k)_[i_1]..._[i_N]
+                       A..._[j_1]..._[j_k])
                   (map 1#g~[i_%1]~[j_%1] (between 1 k))))))))
 
 (define $δ
diff --git a/test/UnitTest.hs b/test/UnitTest.hs
--- a/test/UnitTest.hs
+++ b/test/UnitTest.hs
@@ -31,7 +31,7 @@
         assertEgisonM m = fromEgisonM m >>= assertString . either show (const "")
     
         collectDefsAndTests (Define name expr) (bindings, tests) =
-          ((show name, expr) : bindings, tests)
+          (((stringToVar $ show name), expr) : bindings, tests)
         collectDefsAndTests (Test expr) (bindings, tests) =
           (bindings, expr : tests)
         collectDefsAndTests _ r = r
