Freitag, 18. April 2014

Optimal left to right pattern matching Automata in literate clojure

Optimal left to right pattern matching Automata
Hi, I am in the process of rewriting expresso's rule engine. Currently, rules in expresso can be very expressive, but the rule engine was written in the first two weeks of my gsoc project time on top of core.logic, so there are many efficiency issues with it. I did some research on how the rule engines from other term rewriting engines are build (like maude, elan, stragego ...). Following a few references, I came to this paper, which presents an algorithm to match an expression to multiple simultaneous patterns without backtracking in one pass over the expression, which is really cool and the perfect basis for a good rule engine/compiler. I implemented the algorithm in the paper and also build a rule compiler, that unrolls the automaton that is constructed into an efficient clojure function. It is a literate program, you can find it here. The rest of this post is the created html export from the literate org file.

1 Introduction

This is a clojure implementation of the algorithm from this paper. The problem that this addresses is matching a term against multiple patterns simultaneously, without backtracking scanning the expression only one time from left to right. The patterns are assumed to be linear, that is that there are no two variables with the same name in a pattern. Testing for equality has to be done as an additional check when the linear pattern matched.
To accomplish this, a directed acyclic graph representing a deterministic automaton is created from the pattern, with transitions labeled by the next symbol read during a left to right scan through the pattern and the currently matching patterns as nodes. The dag is kept minimal, that is there are no two states in the dag that produce the same matching sub-tree.
I extended the algorithm in the paper to also work when there is a wildcard on a function symbol like in the following pattern: '(? a b) and also to handle functions with multiple arities. This adds a few extra cases to the interpreter and the compiler, but in the case it isn't needed doesn't slow down the matching process.
Interpreting it works as expected - scan through the input expression, for each symbol follow the labeled transition if it exists - pick the default route if one exists in case that fails - fail otherwise - repeat until at failure state or the end of the expression is reached
The dag can also be compiled to an optimized clojure function resembling the decision tree that the dag represents. Basically, the function consists of a bunch of (case <location-in-expression> <symbol1> <forward-location-and-check-next-input> …. <default> <go-through-default-route-if-possible>) thus eliminating the need to search through the dag at matching time.

1.1 Implementation

(ns optimal-left-to-right-pattern-matching-automata.core
  (:require [clojure.set :as set]
            [clojure.walk :as walk]
            [clojure.zip :as zip]))
We need a (meta-) symbol for a default transition. It will be called omega from now on
(def omega '?)

1.1.1 Representing patterns

Because we are concerned with scanning expressions from left to right, the matching positions of the patterns can be totally ordered - by how right they appear in the printed representation - and put in a single list. Function symbols are represented as [<function-symbol> <number-of-arguments>], so that the flat representation retains all information about the original structure of the pattern. For example, the pattern '(f (g a b) a a b) can be represented as '([f 4] [g 2] a b a a b). In this representation, a pattern is just a list of transition labels that the automaton must perform in order to match an expression against the pattern. During the matching, there will always be a current state which is all patterns with the same matching prefix, a current position where the next symbol will be read, and a suffix to be matched for the next symbols read. This is the definition of a matching-item.
(defn matching-item
  "A matching item is a triple r:a*b where ab is a term and r is a rule label.
   The label identifies the origin of the term ab and hence, in a term rewriting system, the rewrite rule which has to be applied when ab is matched * is called
  the matching dot, a the prefix and b the suffix. The first symbol of b is the matching symbol. The position of the matching dot is the matching position"
  [r a b]
  [r a b])

(defn matching-symbol [matching-item]
  (let [[r a b] matching-item]
    (first b)))

(def infinity (Double/MAX_VALUE))

(defn final? [matching-item]
  (let [[r a b] matching-item]
    (empty? b)))

(defn matching-position [matching-item]
  (if (final? matching-item)
    infinity
    (let [[r a b] matching-item]
      (inc (count a)))))
(defn initial-matching-item [label pattern]
  [label '() pattern])
The current state of the automaton is then a set of matching items which share the same prefix.
(defn matching-set? [matching-items]
  (let [[r a b] (first matching-items)]
    (every? #(let [[r2 a2 b2] %] (= r r2)) (rest matching-items))))

(defn initial-matching-set? [matching-items]
  (and (matching-set? matching-items) (= '() (second (first matching-items)))))

(defn final-matching-set? [matching-items]
  (and (matching-set? matching-items) (= '() (nth (first matching-items) 2))))

1.1.2 The Transition function for the automaton

After we know what the states and the transitions of the automaton will be, we can start looking at the definition for the transition function delta. For more explanation, see the paper itself. Basically, from the current state - the current-matching-set - it returns as next node the set of matching items which could be forwarded by the symbol s - that is what the accept function does. It also avoids backtracking by adding more states when there is an ambiguity in the form that one pattern has a default next transition and another has a transition that goes a level deeper with a function symbol. If the function symbol transition would be followed, it could be that it failed and one had to backtrack and go through the omega transition. Therefore, for each such situation a new pattern is added to the matching set which consists of the omega rule but with the omega replaced by the next function symbol and a number of omegas that match the functions arity. It is also important to do this closing over the current matching set at the very beginning to handle the case of a default omega pattern. The paper fails to mention that.
(defn forward-matching-position [matching-item]
  (let [[r a b] matching-item]
    [r (concat a [(first b)]) (rest b)]))

(defn functions [matching-set]
  (into #{} (filter #(or (and (symbol? %) (not= omega %))
                         (and (sequential? %)
                              (symbol? (first %))
                              (number? (second %))))
                    (map matching-symbol matching-set))))

(defn arity [function-symbol]
  (or (and (sequential? function-symbol) (second function-symbol)) 0))

(defn accept [matching-items s]
  (map forward-matching-position
       (filter #(= (matching-symbol %) s) matching-items)))

(defn close [matching-items]
  (let [F (functions matching-items)]
    (set/union matching-items
               (for [matching-item matching-items
                     function-symbol F
                     :let [arityf (arity function-symbol)]
                     :when (and (= omega (matching-symbol matching-item)))]
                 (let [[r a b] matching-item]
                   [r a (concat [function-symbol] (repeat arityf omega)
                                (rest b))])))))

(defn delta [matching-items s]
  (close (accept matching-items s)))

1.1.3 Creating the DAG

  1. Graph implementation
    Here is a very simple implementation of a functional graph data structure
    ;;quick and dirty functional graph implementation
    (def empty-graph {})
    
    (defn add-node [g n]
      (if (g n)
        g
        (assoc g n {:next #{} :prev #{}})))
    
    (defn add-edge [g n1 n2 l]
      (-> g
          (add-node n1)
          (add-node n2)
          (update-in [n1 :next] conj [n2 l])
          (update-in [n2 :prev] conj [n1 l])))
    
    (defn remove-edge [g n1 n2 l]
      (-> g
          (add-node n1)
          (add-node n2)
          (update-in [n1 :next] disj [n2 l])
          (update-in [n2 :prev] disj [n1 l])))
    
    (defn remove-node [g n]
      (if-let [{:keys [next prev]} (g n)]
        ((comp
          #(dissoc % n)
          #(reduce (fn [g* [n* l*]] (remove-edge g* n* n l*)) % prev)
          #(reduce (fn [g* [n* l*]] (remove-edge g* n n* l*)) % next))
         g)
        g))
    
  2. Recognizing equivalent states
    To make the created automaton minimal, equivalent states have to be recognized during the construction phase. Two states are equivalent, if for each item in set1 there exists an equivalent item in set2. Two matching items are equivalent, if they have the same rule label and the same suffix.
    (defn equivalent-matching-items? [matching-item1 matching-item2]
      (let [[r1 a1 b1] matching-item1 [r2 a2 b2] matching-item2]
        (and (= r1 r2) (= b1 b2))))
    
    (defn extract-first-by
      "returns [extracted rest-of-collection] or false"
      [f coll]
      (loop [[c & cs] coll rest-coll []]
        (if c
          (if (f c)
            [c (concat rest-coll cs)]
            (recur cs (conj rest-coll c)))
          false)))
    
    (defn equivalent-matching-sets? [matching-set1 matching-set2]
      (loop [[mit & mits] matching-set1 matching-set2 matching-set2]
        (if mit
          (if-let [[mit2 mits2] (extract-first-by #(equivalent-matching-items? mit %)
                                                  matching-set2)]
            (recur mits mits2)
            false)
          (empty? matching-set2))))
    
  3. Constructing the DAG
    For detailed description about this algorithm, see the paper. Basically, we start with the initial-matching-set and create new states for all possible transitions, add the nodes and the edges to the graph, or only the transition if there already exists an equivalent state in the graph. Then sort the newly created states according to their matching position, so that states with only a few already matched items are handled first. The creation ends when the list of states is traversed completely.
    (defn failure? [state]
      (or (= '() state) (nil? state)))
    
    (defn get-next-node [g n l]
      (some #(and (= (second %) l) (first %)) (get-in g [n :next])))
    
    (defn search-equivalent-node [graph node]
      (first (for [[n v] graph
                   :when (equivalent-matching-sets? node n)]
               n)))
    
    (defn insert-according-to-matching-position [nodes-to-visit new-matching-set]
      ;;nodes-to-visit has to be sorted according to matching-position
      ;;all matching positions in a matching set are the same
      (let [nmp (matching-position (first new-matching-set))]
        (loop [[n & ns :as nodes-left] nodes-to-visit new-nodes-to-visit []]
          (if n
            (if (<= (matching-position (first n)) nmp)
              (recur ns (conj new-nodes-to-visit n))
              (vec (concat new-nodes-to-visit [new-matching-set] nodes-left)))
            (conj nodes-to-visit new-matching-set)))))
    
    ;;problem hier? gibt nur ein omega jetzt mehrere
    (defn create-new-states [pos nodes-to-visit graph]
      (let [current-state (nth nodes-to-visit pos)
            F (functions current-state)]
        (loop [[s & ss] (concat F [omega]) nodes-to-visit nodes-to-visit graph graph]
          (if s
            ;;work to do
            (let [new-matching-set (delta current-state s)]
              ;;check if there is already an equivalent matching-set in the graph
              (if-let [eq-node (search-equivalent-node graph new-matching-set)]
                (recur ss nodes-to-visit (add-edge graph current-state eq-node s))
                (recur ss (insert-according-to-matching-position
                           nodes-to-visit new-matching-set)
                       (add-edge graph current-state new-matching-set s))))
            ;;all symbols consumpted, so return the new state
            [graph nodes-to-visit]))))
    
    (defn create-dag [initial-matching-set]
      (loop [graph empty-graph nodes-to-visit [initial-matching-set] pos 0]
        (if (= (count nodes-to-visit) pos)
          ;;all nodes visited, so return graph
          (remove-node graph '())
          (let [[new-graph new-nodes-to-visit]
                (create-new-states pos nodes-to-visit graph)]
            (recur new-graph new-nodes-to-visit (inc pos))))))
    

1.1.4 Interpreting the DAG

With the constucted minimal dag like described in the paper, we can leave it and now implement how to interpret that dag to match an expression against multiple paterns. To do this, we will traverse the expression from left to right using clojure zippers. We recursively check for the next transition, follow it and move the zipper forward accordingly and fail if there is no transition possible. If we go through a wildcard then we add the current value of the zipper location to the bindings ;;TODO may miss some bindings in rules created by close
(defn consume-next [g current-state symbol]
  (let [next-state (get-next-node g current-state symbol)]
    (if (failure? next-state)
      ;;there was no link, so go through omega link
      [(get-next-node g current-state omega) [symbol]]
      [next-state []])))

(defn consume-next-level-down [g current-state [symbol count]]
  (let [next-state (get-next-node g current-state [symbol count])]
    (if (failure? next-state)
      ;;there was no link, so go through omega link
      [(get-next-node g current-state [omega count]) [symbol]]
      [next-state []])))


(defn- next-without-down
  [loc]
  (if (= :end (loc 1))
    loc
    (or
     (zip/right loc)
     (loop [p loc]
       (if (zip/up p)
         (or (zip/right (zip/up p)) (recur (zip/up p)))
         [(zip/node p) :end])))))

(defn match-expression [g patterns expression]
  (loop [loc (zip/seq-zip expression) node patterns bindings []]
    (if (or (failure? node) (zip/end? loc))
      ;;done
      [node bindings]
      (if (zip/branch? loc)
        ;;ok try if head symbol matches
        ;;we are using preorder throughout matching
        (let [children-count (dec (count (zip/children loc)))
              head-loc (zip/next loc)
              [next-node add-bindings]
              (consume-next-level-down g node [(first head-loc) children-count])]
          (if (failure? next-node)
            ;;head got no match so we have to stay at the original level and try
            ;;to match there for a value or omega
            (let [[next-node add-bindings] (consume-next g node (first loc))]   
              (recur (next-without-down loc) next-node
                     (concat bindings add-bindings)))
            ;;head location got a match so we go on on this level
            (recur (zip/next head-loc) next-node
                   (concat bindings add-bindings))))
        ;;we have no possibility to go down a level deeper so we can just
        ;;consume directly
        (let [[next-node add-bindings] (consume-next g node (first loc))]
          (recur (zip/next loc) next-node
                 (concat bindings add-bindings)))))))
  1. Testing
    Here are a few sample calls and tests:
    (use 'clojure.test)
    (let [initial-matching-set (close [(initial-matching-item 1 '([? 2] a b))
                                       (initial-matching-item 2 '([? 1] a))
                                       (initial-matching-item 3 '(?))])
          dag (create-dag initial-matching-set)]
      (is (= '[([3 (?) ()]) (1)] (match-expression dag initial-matching-set 1)))
      (is (= '[([3 ([? 1] a) ()] [2 ([? 1] a) ()]) (+)]
             (match-expression dag initial-matching-set '(+ a))))
      (is (= '[([3 ([? 2] a b) ()] [1 ([? 2] a b) ()]) (+)]
             (match-expression dag initial-matching-set '(+ a b))))
      (is (= '[([3 (?) ()]) ((+ a b c))]
             (match-expression dag initial-matching-set '(+ a b c)))))
    

1.1.5 Compiling the DAG to a fast clojure function

The expression matching can be taken a level further, to the point that the dag can be compiled to a fast clojure function. The resulting clojure function will look like this:
(fn [expression]
       (let [loc (zip/seq-zip expression)]
         ;;now code for the single transitions
         (or 
           ;;if there are possible transitions in the dag that lead one
           ;;level down - if now the next part is replaced by false
           ;;and the next branch of the or is taken
           (and (zip/branch? ~'loc) ;;fail if we are not in a branch
              (let [head-loc (zip/next loc)]
                (case (first head-loc) ;;fast dispatch on the function symbol
                  <function-symbol> (and (check-if-argument-count-matches)
                                         <code-for-the-next-transitions>)
                  #_ (...)
                  ;;default case
                  <code-for-wildcard-transition or nil if no wildcard>)))
           ;;if there is no matching transition for the current head symbol
           ;;try matching the whole current subtree
           (case (first loc)
             <variable-or-constant> <code-for-next-transition>
             #_ (...)
             <code-for-wildcard-transition or nil if there is no wildcard>))))
In the end-nodes of the decision tree the code returns either nil for a failure node or sorts the applicable rules by priority (currently only their label but one could introduce the rule that more specific rules come first) and for each defines the bindings, checks the conditions and returns their result.
Therefor, we now extend the notion of a pattern to the notion of a rule. Currently this is really low level and the rule engine on top if this should take a more human readable form.
A rule has the form [<label> <pattern> <conditions> <results> <wildcard-positions>] label and pattern are the same as before, conditions is just a list of expressions to evaluate after succesful match, result is the rhs of the rule and wildcard-positions maps the wildcards in the pattern to the positions in the expression.
With this the compile-rules function can be defined
(defn get-in-expression [expression key]
  (loop [loc (zip/seq-zip expression) [k & ks] key]
    (if k
      (let [new-loc (loop [k k loc (zip/down loc)]
                      (if (> k 0)
                        (recur (dec k) (zip/right loc))
                        loc))]
        (recur new-loc ks))
      (first loc))))

(defn compile-step [g current-state rule-map]
  (let [possible-moves (doall (map last (:next (get g current-state))))
        head-moves (doall (filter sequential? possible-moves))
        current-level-moves (doall (remove sequential? possible-moves))]
    (if (empty? possible-moves) 
      `(and (zip/end? ~'loc)
            ;;current-state was successfully matched. Now get the results for the
            ;;matched rules in current-stater
            (or
             ~@(for [[label & rest] (sort-by first (filter final? current-state))
                     :let [[conditions result omga-positions]
                           (get rule-map label)]]
                 `(let [~@(for [[name pos] omga-positions
                                entry
                                [name `(get-in-expression ~'expression ~pos)]]
                            entry)]
                    (and ~@(concat conditions [result]))))))
      `(or ~(if (empty? head-moves)
              'false
              ;;have to test going a level deeper
              `(and (zip/branch? ~'loc)
                    (let [~'head-loc (zip/next ~'loc)]
                      (case (first ~'head-loc)
                        ;;now all next steps have to be written down in a
                        ;;case - the right hand side will be a recursive
                        ;;call to create the code at the next level
                        ;;the default of case is either nil or the level
                        ;;from following a [? <number>] label in the graph
                        ~@(concat
                           (for [[s c] head-moves :when (not= omega s)
                                 entry
                                 [s `(and
                                      (= (dec (count (zip/children ~'loc))) ~c)
                                      (let [~'loc (zip/next ~'head-loc)]
                                        ~(compile-step
                                          g
                                          (get-next-node
                                           g current-state [s c])
                                          rule-map)))]]
                             entry)
                           [(let [omega-downs (filter #(= (first %) omega)
                                                      head-moves)]
                              `(case (dec (count (zip/children ~'loc)))
                                 ~@(concat
                                    (for [[omga c] omega-downs
                                          entry
                                          [c
                                           `(let [~'loc (zip/next ~'head-loc)]
                                              ~(compile-step
                                                g
                                                (get-next-node
                                                 g current-state[omega c])
                                                rule-map))]]
                                      entry)
                                    ;;no further defaulting possible - fail
                                    '(nil))))])))))
           (case (first ~'loc)
             ~@(concat
                (for [symbol current-level-moves :when (not= omega symbol)
                      entry
                      [symbol `(let [~'loc (next-without-down ~'loc)]
                                 ~(compile-step
                                   g
                                   (get-next-node g current-state symbol)
                                   rule-map))]]
                  entry)
                [(if (some #{omega} current-level-moves)
                   ;;we have a default case to fall back to
                   `(let [~'loc (next-without-down ~'loc)]
                      ~(compile-step
                        g
                        (get-next-node g current-state omega)
                        rule-map))
                   'nil)]))))))



(defn compile-rules [rules]
  (let [res
        (for [[label pattern conditions result omga-positions] rules]
          [(initial-matching-item label pattern) [label [conditions result
                                                         omga-positions]]])
        initial-matching-set (close (map first res))
        rule-map (into {} (map second res))
        dag (create-dag initial-matching-set)]
    `(fn [~'expression]
       (let [~'loc (zip/seq-zip ~'expression)]
         ~(compile-step dag initial-matching-set rule-map)))))
  1. Tests with example rules
    Here are two example rules: (f a a ?a a) => ?a (f (g a ?b) a ?b a) => ?b Encoded in the current low-level representation they become
    [[1 '([f 4] a a ? a) [] '?a '{?a [3]}]
      [2 '([f 4] [g 2] a ? a ? a) '[(= ?a ?b)] '?b '{?b [1 2] ?a [3]}]]
    
    Here are the corresponding tests:
    (let [rules
          [[1 '([f 4] a a ? a) [] '?a '{?a [3]}]
           [2 '([f 4] [g 2] a ? a ? a) '[(= ?a ?b)] '?b '{?b [1 2] ?a [3]}]]
          f (eval (compile-rules rules))]
      (is (= 'c (f '(f (g a c) a c a))))
      (is (not (f '(f (g a b) a c a))))
      (is (= 'a (f '(f a a a a))))
      (is (not (f '(f a a a b)))))
    
  2. Example code
    The compiled code for the two rules above looks like this:
    (clojure.core/fn
     [expression]
     (clojure.core/let
      [loc (clojure.zip/seq-zip expression)]
      (clojure.core/or
       (clojure.core/and
        (clojure.zip/branch? loc)
        (clojure.core/let
         [head-loc (clojure.zip/next loc)]
         (clojure.core/case
          (clojure.core/first head-loc)
          f
          (clojure.core/and
           (clojure.core/=
            (clojure.core/dec
             (clojure.core/count (clojure.zip/children loc)))
            4)
           (clojure.core/let
            [loc (clojure.zip/next head-loc)]
            (clojure.core/or
             (clojure.core/and
              (clojure.zip/branch? loc)
              (clojure.core/let
               [head-loc (clojure.zip/next loc)]
               (clojure.core/case
                (clojure.core/first head-loc)
                g
                (clojure.core/and
                 (clojure.core/=
                  (clojure.core/dec
                   (clojure.core/count (clojure.zip/children loc)))
                  2)
                 (clojure.core/let
                  [loc (clojure.zip/next head-loc)]
                  (clojure.core/or
                   false
                   (clojure.core/case
                    (clojure.core/first loc)
                    a
                    (clojure.core/let
                     [loc
                      (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                       loc)]
                     (clojure.core/or
                      false
                      (clojure.core/case
                       (clojure.core/first loc)
                       (clojure.core/let
                        [loc
                         (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                          loc)]
                        (clojure.core/or
                         false
                         (clojure.core/case
                          (clojure.core/first loc)
                          a
                          (clojure.core/let
                           [loc
                            (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                             loc)]
                           (clojure.core/or
                            false
                            (clojure.core/case
                             (clojure.core/first loc)
                             (clojure.core/let
                              [loc
                               (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                                loc)]
                              (clojure.core/or
                               false
                               (clojure.core/case
                                (clojure.core/first loc)
                                a
                                (clojure.core/let
                                 [loc
                                  (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                                   loc)]
                                 (clojure.core/and
                                  (clojure.zip/end? loc)
                                  (clojure.core/or
                                   (clojure.core/let
                                    [?b
                                     (optimal-left-to-right-pattern-matching-automata.core/get-in-expression
                                      expression
                                      [1 2])
                                     ?a
                                     (optimal-left-to-right-pattern-matching-automata.core/get-in-expression
                                      expression
                                      [3])]
                                    (clojure.core/and
                                     (= ?a ?b)
                                     [?a ?b])))))
                                nil))))))
                          nil))))))
                    nil))))
                (clojure.core/case
                 (clojure.core/dec
                  (clojure.core/count (clojure.zip/children loc)))
                 nil))))
             (clojure.core/case
              (clojure.core/first loc)
              a
              (clojure.core/let
               [loc
                (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                 loc)]
               (clojure.core/or
                false
                (clojure.core/case
                 (clojure.core/first loc)
                 a
                 (clojure.core/let
                  [loc
                   (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                    loc)]
                  (clojure.core/or
                   false
                   (clojure.core/case
                    (clojure.core/first loc)
                    (clojure.core/let
                     [loc
                      (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                       loc)]
                     (clojure.core/or
                      false
                      (clojure.core/case
                       (clojure.core/first loc)
                       a
                       (clojure.core/let
                        [loc
                         (optimal-left-to-right-pattern-matching-automata.core/next-without-down
                          loc)]
                        (clojure.core/and
                         (clojure.zip/end? loc)
                         (clojure.core/or
                          (clojure.core/let
                           [?a
                            (optimal-left-to-right-pattern-matching-automata.core/get-in-expression
                             expression
                             [3])]
                           (clojure.core/and ?a)))))
                       nil))))))
                 nil)))
              nil))))
          (clojure.core/case
           (clojure.core/dec
            (clojure.core/count (clojure.zip/children loc)))
           nil))))
       (clojure.core/case (clojure.core/first loc) nil))))
    
Author: Maik Schünemann
Created: 2014-04-18 Fri 12:51
Emacs 24.2.1 (Org mode 8.2.5g)

Montag, 23. September 2013

Finished GSoC project Expresso

Finished GSoC project Expresso

Finished GSoC project Expresso


GSoC ends today and I can announce the 0.2.0 version of the expresso library. It is build on top of core.logic and core.matrix and provides symbolic manipulation of algebraic expressions.

What's there?

  1. An api/dsl for manipulation of algebraic expressions which doesn't get in your way. Expresso's expressions are just clojure s-expressions and can be manipulated with rich set of clojure sequence functions
  2. useful manipulations for mathematical expressions: simplify, multiply-out, differentiate, …
  3. An equation solver which is capable of solving a single equation and multiple equations for unknowns.
  4. An optimizer which transforms a mathematical expression to a semantically equivalent but performanter one
  5. An expression compiler to compile an expression to an efficient clojure function
  6. A semantic rule based translator on top of which many of expresso's features are implemented

The code is fully documented and I wrote a tutorial and showcase of expresso, the expresso-tutorial.

GSoC has been a really fun and valuable time for me. I learned a lot. Of course I will continue developing expresso! Expresso and core.matrix are the first steps in the direction of a full computer algebra system for clojure. I hope that it will help clojure to be an attractive choice for scientific computing projects in the future.

Showcase: Here are two examples of expresso's facility to manipulate mathematical expressions. They can be found and are explained in the expresso-tutorial 1.

  1. solving word problems:
(solve 'blue
  (ex (= pencils (+ green white blue red)))
  (ex (= (/ pencils 10) green))
  (ex (= (/ pencils 2) white))
  (ex (= (/ pencils 4) blue))
  (ex (= red 45))) ;=> #{{blue 75N}}
  1. Analyzing roots and extremata of functions. This code shows how easy one can implement tasks involving symbolic manipulation with expresso:
(defn roots
  "returns the set of roots of the expression in regard to var"
  [var expr]
  (solve var (ex (= ~expr 0))))


(defn extremata 
  "gets the extrema of the expression in regard to var. Returns a map with the
   keys :maxima and :minima"
  [var expr]
  (let [d1 (differentiate [var] expr)
 d2 (differentiate [var] d1)
 candidates (roots var d1)]
    (if (seq candidates)
      (let [extremata
     (->> candidates
   (map (fn [candidate] [candidate (evaluate d2 {var candidate})]))
   (remove #(== 0 (second %)))
   (group-by #(< 0 (second %))))]
 {:maxima (map first (get extremata false))
  :minima (map first (get extremata true))}))))


(defn analyse-function 
  "returns a map with the :roots, the :maxima and the :minima of the expression
   in regard to var"
  [var expr]
  (assoc (extremata var expr)
    :roots (roots var expr)))

(analyse-function 'x (ex (- (** x 4) (** x 2))))
;=> {:roots #{0 -1 1},
;;   :maxima (0),
;;   :minima (0.7071067811865476 -0.7071067811865476)}

Ideas/feedbacks etc are greatly appreciated! Enjoy, Maik

Montag, 26. August 2013

first release of expresso

First release of expresso

First release of expresso


It is now ca 1 month before gsoc ends. It has been a fantastic time and valuable learning experience so far!

Expressos basic functionality is now in place. It is convenient for handling expressions, can simplify/differentiate/solve/optimize them and compile them to efficient clojure functions.

Maybe there are some early adopters, who want to use expresso's symbolic manipulation and want to give it a try! (and report issues/bugs etc) See the github repository and README.md for details.

Samstag, 27. Juli 2013

symbolic expression manipulation with expresso

symbolic expression manipulation with expresso

symbolic expression manipulation with expresso


It is now half way through the google summer of code project and I will do a quick showoff in this post about what is supported in expresso right now and what will come in the near future.

1 Constructing expressions

Expresso expressions are clojure s-expressions. Even if expresso can use custom types to better represent for example poynomials internally they will implement ISeq so that they can be used like regular s-expressions in every regard. Expresso provides a few convenient methods to help you construct expressions:

(use 'numeric.expresso.core)

;; the ex macro is convenient for quickly constructing expressions. variables 
;; are automatically quoted and the function symbols are automatically namespace 
;; qualified, which is important because clojure.core/* behaves different than 
;; core.matrix.operators/*

(ex (+ x y (* x 2 z))) ;=> (clojure.core/+ x y (clojure.core/* x 2 z))

;;If you have a local named x and would like to have the value of x in the 
;;expression you need to unquote it
(let [x 3]
  (ex (+ 1 ~x))) ;=> (clojure.core/+ 1 3)

;;If you have more symbols for which you want the actual values (for example if
;; the constants aren't primitive) Ther is the ex' macro. You have to explicitly 
;; ' (quote) the arguments with it therefore the name :)
(let [x 3]
  (ex' (+ x 'x))) ;=> (clojure.core/+ 3 x)
;; It also takes an optional symbol vector as first argument. The variables in it
;; will be implicitly quoted
(ex' [x] (* 3 x)) ;=> (clojure.core/* 3 x)

;; these macros compile down to calls the the more fundamental ce function 
;; (construct expression) which
;; creates an expression from the symbol and arguments.
(ce `+ [1 2 3]) ;=> (clojure.core/+ [1 2 3])
;;you can think of it like list* + adding information expresso uses in the 
;;algebraic manipulation.

;;construct-with lets you transform your normal heavy number crunshing code to
;; clean symbolic expression generating code. In its body it replaces all the
;; symbols contained in the symbol vector with the equivalent expresso construcing 
;;functions. Is is especially convenient for writing rules.
(construct-with [+ - * /]
  (defn create-silly-expression [x]
    (+ 1 (* 2 (- 3 (/ 4 x))))))

(create-silly-expression 5)
;=>(clojure.core/+ 1 (clojure.core/* 2 (clojure.core/- 3 (clojure.core// 4 5))))

2 Basic manipulations of expressions

You can evaluate an expression providing a symbol map for the symbols in it and you can substitute partsof the expression

(evaluate (ex (+ (* 2 x) (* x 2))) {'x 2}) ;=> 8
(eval (substitute (ex (+ (* 2 x) (* x 2))) {'x 2})) ;=> 8

(substitute (ex (+ (* 1 0) 4)) {(ex (* 1 0)) 0}) ;=> (clojure.core/+ 0 4)

3 A semantic rule based translator on top of core.logic

Term rewriting through rule application is a - or 'the' - fundamental technique for all sorts of expression manipulation. Therfore, expresso would not live up to its aim to be a general clojure library for manipulating expressions without a powerful rule based translator. Here is what the rule based translator looks like:

(use 'numeric.expresso.core)

;;rules are constructed with the rule macro. Rules in expresso are semantic by 
;;default, that means commutative expression match regardlessy of their order of 
;;arguments for example. The syntax is (rule pat :=> trans) where pat is an 
;;expression which can contain variables and trans can be an expression or a core.logic
;;relation which will be called upon succesful matching. An optional guard relation can
;;be specified at the end whith :if guard. See my previous posts for more examples

(construct-with [+ - * /]

(def rules [(rule (+ 0 ?x) :=> ?x)
            (rule (+ 0 ?&*) :=> (+ ?&*))])

(apply-rule (first rules) (+ 0 1)) ;=> 1
(apply-rule (first rules) (+ 1 0)) ;=> 1
(apply-rule (first rules) (+ 0 1 2)) ;=> nil
(apply-rule (second rules) (+ 0 1 2)) ;=>(clojure.core/+ 1 2)
(apply-rule (second rules) (+ 1 0 3 4 2)) ;=> (clojure.core/+ 1 3 4 2)

;;apply-rule tries to apply the given rule. apply-rules tries to apply a whole 
;;vector of rules and transform-expression transforms a given expression according
;;to a rule vector to a normal form on which no rules can be applied anymore.
;;the ?&* represents a sequence-matcher, which matches a subexpression and spits 
;;it in after matching. This allows powerful transformations to be represented as rules
;;advanced rules:
;;The rule macro also supports to define a transformation function inline,
;; with the ;==> syntax and extractors, which can take apart an expression 
;;using a custom extracting function. Together they enable really powerful 
;;translations. For example the following could be part of a function which 
;;rearranges an expression in regard to a variable v using the rule based 
;;translator like this
(with-expresso [+ cons? = - / * ]
(def rearrange-rules
  [(rule [(cons? ?p ?ps) (= (+ ?&+) ?rhs)]
         :==> (let [[left x right] (split-in-pos-sm ?&+ ?p)]
                [?ps (= x (- ?rhs left right))]))
   (rule [(cons? ?p ?ps) (= (* ?&+) ?rhs)]
         :==> (let [[left x right] (split-in-pos-sm ?&+ ?p)]
                [?ps (= x (/ ?rhs (* left right)))]))]))

;;cons? is an extractor which matches a vector and binds ?p to the first and ?ps to
;; the rest of it the part after the :==> is normal cojure code!
;;it is also possible to define custom extractors

4 algebraic expression manipulation algorithms

This is the second big part of expresso which contains functions to solve and rearrange an expression, to simplify it etc. The range of expressions these functions can handle will be greatly enhanced in the second part of gsoc

(simplify (ex (+ (* x 0) (* x (- x x))))) ;=> 0
(simplify (ex (+ (* 2 x) 4 5 (* 3 x)))) ;=> (clojure.core/+ (clojure.core/* 5 x) 9)

;;differentiate takes the derivative of the given expression regarding the
;;variable v
(differentiate 'x (ex (+ (* 3 (** x 2)) (* -2 x))))
;=>(clojure.core/+ (clojure.core/* 6 x) -2)
(differentiate 'y (ex (+ (* 3 (** x 2)) (* -2 x))))
;=> 0

;;rearrange rearranges an expression containing one occurrence of the specified 
;;variable to the form (= v rhs)
(rearrange 'x (ex (= (+ 1 x) 4))) ;=> (clojure.core/= x (clojure.core/- 4 1))
(rearrange 'x (ex (= (+ 1 (* 2 (+ 3 (* 4 x)))) 10)))
;=> (clojure.core/= x (clojure.core// (clojure.core/- 
;   (clojure.core// (clojure.core/- 10 1) (clojure.core/* 2)) 3) (clojure.core/* 4)))

;;solve tries to simplify the expression so that it contains only one occurrence of 
;;the variable, rearranges it then and simplifies the rhs
(solve 'x (ex (= (+ 1 x) 3))) ;=> (clojure.core/= x 2)
(solve 'x (ex (= x (+ x 1)))) ;=> () ;no solution found
(solve 'x (ex (= x x))) ;=> _0 ; x undetermined
(solve 'x (ex (= (* 3 x) (+ (* 4 x)  4)))) ;=> (clojure.core/= x -4)

5 optimizing expressions

Like the solving part above this is in focus for the second part of gsoc.

;;optimize takes an expression optimizes it and returns a function which can be 
;;called with the symbols in the expression.

(def compiled-expr (optimize (+ (* 1 x) (* x 1))))

(compiled-expr {'x 2}) ;=> 4

;;optimize doesn't do very much by now but it recognises common subexpressions
;; and only executes them once to see what optimize has done lets check the result
;; of optimize* which doesn't return a callable function but the optimized expression


(def optimized (optimize* (ex (+ (* 1 x) (* x 1)))))

optimized ;=> (let [local50510 (clojure.core/* 1 x)] 
               ;(clojure.core/+ local50510 local50510))

(evaluate optimized {'x 2}) ;=> 4

all the code you see here is working in the current master branch of expresso. As always comments/critique/thoughts welcome.

Date: 2013-07-29T14:43+0200

Author: Maik Schünemann

Org version 7.8.09 with Emacs version 23

Validate XHTML 1.0

Samstag, 29. Juni 2013

Sequential Matching in Expresso

Sequential Matching in Expresso

Sequential Matching in Expresso


I In the second week of gsoc, I extended the rule syntax to handle sequential matchers, to have proper support for functions having variadic arguments, what is quite common in clojure. They, as the name suggest, match zero or more arguments in the expression. In this post, I like to demonstrate their use and how they enable rules to be more expressive.

1 Introduction

For introduction, here is an expresso rule that cancels out a zero in an expression

(use 'numeric.expresso.rules)
(use 'numeric.expresso.construct)
(with-expresso [+]

  (def cancel-zero-rule (rule (+ 0 ?x) :=> ?x))

  (apply-rule cancel-zero-rule (+ 3 0)) ;=> 3
  (apply-rule cancel-zero-rule (+ 3 0 2)) ;=> nil
)

Because expresso knows that clojure.core/+ is a commutative operation, it matches not regarding the order of arguments, but if there is a permutation, which will match. Therefore, the first call to apply-rule succeeds binding ?x to 3. But, as the rule is written, it expects + to take exactly two arguments. That is why the secand call to apply-rule fails. In many languages, + takes two arguments, but in clojure (and other lisps) + (and many other functions) are variadic, so they support zero, one or more arguments. Below is the version of the rule, which supports the lispy variadic + function

(use 'numeric.expresso.rules)
(use 'numeric.expresso.construct)
(with-expresso [+]

  (def cancel-zero-rule (rule (+ 0 ?&*) :=> (+ ?&*)))

  (apply-rule cancel-zero-rule (+ 3 0)) ;=> (clojure.core/+ 3)
  (apply-rule cancel-zero-rule (+ 3 0 2)) ;=>(clojure.core/+ 3 2) 
)  

The commutative matching functions only allows for one seq-matcher in the pattern. More than one seq-matcher would not make sense, if the order of the arguments doesn't matter. The normal expression match function supports multiple seq-matchers, I will give an example for that below.

2 Example: Simplification Rules in Expresso

To demonstrate the usefulness of expressos rule syntax I present a few common simplification rules for + and *

(with-expresso [+ - * /]

  (def simplification-rules
    [(rule (+ 0 ?&*) :=> (+ ?&*))
     (rule (+ ?x ?x ?&*) :=> (+ (* 2 ?x) ?&*))
     (rule (* 0 ?&*) :=> 0)
     (rule (* 1 ?&*) :=> (* ?&*))
     (rule (+ (* ?x ?&*a) (* ?x ?&*b) ?&*r)
           :=> (+ (* ?x (+ (* ?&*a) (* ?&*b))) ?&*r))
     (rule (*) :=> 1)
     (rule (+) :=> 0)
     (rule (* ?x) :=> ?x)
     (rule (+ ?x) :=> ?x)])

  (transform-with-rules simplification-rules (+ 1 2 0)) ;=> (clojure.core/+ 1 2)
  (transform-with-rules simplification-rules (* (+ 'x 'x) 1))
   ;=> (clojure.core/* 2 x)
  (transform-with-rules simplification-rules (+ 'x (* 'x 3 1) 'x (* 'y 0)))
   ;=>(clojure.core/* x (clojure.core/+ 3 2))

Especially the fifth rule shows the power of a custom matching function and seq-matchers. It quite literally states the simplification of factoring out. Also note the last call to transform-with-rules. It simplifies (* 'y 0) => 0 and (* 'x 3 1) => (* 'x 3) and (+ 'x (* 'x 3) 'x 0) => (+ 'x (* 'x 3) 'x) => (+ (* 2 'x) (* 'x 3)) => (* 'x (+ 3 2)).

3 Cool example of sequential matching

I want to finish with a very cool example of sequential matching in an associative context. Here, ° denotes a list concatenation operator, so that (° 1 2 3) is the list containing 1 2 3. It is possible to expresso insertion sort as application of only one rule. Here is how:

(use 'clojure.core.logic)
(defn biggero [x y] (project [x y] (== true (> x y))))
(with-expresso [°]

  (def sort-rule (rule (° ?&*1 ?x ?&*2 ?y ?&*3) :=> (° ?&*1 ?y ?&*2 ?x ?&*3)
                       :if (biggero ?y ?x)))

  (transform-with-rules [sort-rule] (° 1 4 2 6 5 4 3 7 8 9))
   ;=> (numeric.expresso.construct/° 9 8 7 6 5 4 4 3 2 1)
)

Btw: I don't recommend using this from a performance perspective :)

4 Whats next

In the next week I want to give the rule based translator a test by writing rules to simplify an expression to a normal form, which could also be the default input for the more sophisticated transformations of expresso. I also want to add support for a ?&+ seq-matcher, which matches one or more arguments

Sonntag, 23. Juni 2013

[GSOC] Expressions, Rules and first Transformations

[GSOC] Expressions, Rules and first Transformations

[GSOC] Expressions, Rules and first Transformations

This is the status report after the 1st week of gsoc. While I still have normal schedule of university classes, I managed to get some time for coding (could'nt stop me doing this)

1 Expressions

Expressions are encoded as s-exp with metadata. Therefore, they can be constructed by hand, or with the expresso.construct namespace for convenience. The central function is ex, which generates the s-exp with op as symbol and args as arguments and with the appropriate metadata for the operation.

(use 'numeric.expresso.construct)

(ex `* 1 2 3) ;=> (clojure.core/* 1 2 3)
(properties (ex `/ 3 2 1)) 
 ;=> [:n-ary [:inverse-of clojure.core/*]]

ex gets the properties of the expression from the multimethod props, which can be extended by users

(defmethod props 'exp [_] [:positive [:arity 1]])

(properties (ex 'exp 0)) ;=> [:positive [:arity 1]]

constructing expressions with ex is still tedious, so there is the with-expression macro.

(with-expresso [+ - * /] (+ 1 2 (* 3 4))) 
;=> (clojure.core/+ 1 2 (clojure.core/* 3 4))
(properties
  (with-expresso [+ - * /] (+ 1 2 (* 3 4)))) 
;=> [:associative :commutative :n-ary]

It takes a vector as first argument and replaces in the code all occurences of the symbols with their namespaced-qualified name and calls ex with it.

Because expressions are just s-expressions with metadata (which gets preserved by core.logic) it is easy to manipulate them with core.logic.

2 Rules

Rules are the most important part of the core expression rule based transformer (obviously), so I really want a good rule design, in which

  • encoding arbitrary transformations as rules is supported
  • it is easy to construct rules syntactically
  • rules are not a black-box for expresso, so that rules can be optimized (i.e. pruning the rules to those applicable to the expression)

2.1 Rule representation

Expresso rules are represented by a vector of [pat trans guard] where pat is an syntactical pattern, which gets matched against the expression. Guard is a core.logic relation, which gets checked after a succesfull match. Trans specifies the new form of expression. In the simple case it can be a syntactic expression, or it can be a core.logic relation (details later). The macro rule constructs a rule. Variables are encoded as a symbol starting with ?.

(use 'numeric.expresso.rules)
(use 'clojure.core.logic)
(with-expresso [* + - = /]
(def rules [(rule (* ?x 1) :=> ?x)
            (rule (* ?x 0) :=> 0)
            (rule (+ ?x 0) :=> ?x)
            (rule (+ ?x (- ?x)) :=> 0)
            (rule (- ?x ?x) :=> (- (* 2 ?x)))])
;if no guard is provided rule inserts the succeed core.logic goal.

(def with-guard (rule (= (* ?a 'x) ?y) :=> (= 'x (/ ?y ?a)) :if (!= ?a 0)))

(apply-rule with-guard (= (* 'x 1) 2)) ;=> (clojure.core/= x (clojure.core// 2 1))
(apply-rule with-guard (= (* 'x 0) 2)) ;=> nil

;the transformation can be a relation itself,
(use 'numeric.expresso.solve)
(def calc-rule (rule ?x :=> (fn [res] (resulto ?x res)) :if (no-variablo ?x)))

(apply-rule calc-rule (* 3 (+ 2 1))) ;=>  9
(apply-rule calc-rule (* 'x 3)) ;=> nil
)

2.2 Expressions specify the matching algorithm

Look closely at the example of the rule with-guard. applying it to (= (* 'x 1) 2) matches with (= (* ?a 'x) ?y). This is not a strange bug, but intented behaviour. Let me explain the reasons for this:

If you have a rule based translator based on normal unification or pattern matching, it will match syntactically, not semantically. This is a problem, because, for example, (* 3 4) and (* 4 3) are semantically equivalent, but not syntactically. So to simplify multiplication with zero you would have to write multiple rules

(with-expresso [*]
  (rule (* ?x 0) :=> 0)
  (rule (* 0 ?x) :=> 0))

It is getting (a lot) worse if you have expressions like (* 2 3 0 1), which is normal in clojure.

Therefore, in expresso, the rules decide the matching algorithms. And again, the matching function is stored in a multimethod (which defaults to syntactical matching) and therefore extensible to new operators.

2.3 Write one, match many

Another situation likely to aries when writing rules for (say) simplification is that you have to duplicate rules with an other operator which behaves the same. Basic example: If you use different implementations for the same function - like clojure.core/+ and an optimized other implementation . you would have to dublicate the rules concerning clojure.core/+ with the other implementation.

Expresso uses clojures ad-hoc hierarchies for this. It matches a pattern against an expression, if the operator of expression is an instance of the operator of the pattern. This enables you to do the following

(derive 'clojure.core/+ e/+)
(derive 'other.ns/+ e/+)

(with-expresso [+ other.ns/+ e/+]
  (def r (rule (e/+ ?x 0) :=> ?x))

  (apply-rule r (+ 3 0)) ;=> 3
  (apply-rule r (other.ns/+ 3 0)) ;=> 3
)

2.4 Going absolutely crazy

This makes it also possible to have very powerful abstract rules. Here is an example for representing mathematical properties of groups. It also shows the function transform-with-rules which recursively applies a vector of rules to the expression


(defmulti group-id identity)
(defmethod group-id :default [_] false)
(defmethod group-id [`+ 0] [_] true)
(defmethod group-id [`* 1] [_] true)

(defn is-group-identityo [op id]
  (project [op id]
           (== true (group-id [op id]))))

(defmulti group-ne identity)
(defmethod group-ne :default [_] nil)
(defmethod group-ne `+ [_] 0)
(defmethod group-ne `* [_] 1)

(defn group-neo [op res]
  (project [op]
           (== res (group-ne op))))

(with-expresso [- /]
(defmulti group-inverse first)
(defmethod group-inverse :default [_] nil)  
(defmethod group-inverse `+ [[_ x]] (- x))
(defmethod group-inverse `* [[_ x]] (/ 1 x))
)
(use 'numeric.expresso.utils)
(defn is-group-inverso [op x y]
  (project [op x]
           (if-let [inv (group-inverse [op x])]
             (== y inv)
             fail)))
(def group-rules
  [(rule (ex ?op ?x ?id) :=> ?x :if (is-group-identityo ?op ?id))
   (rule (ex ?op ?x ?y) :=> (fn [res] (group-neo ?op res)) 
    :if (is-group-inverso ?op ?x ?y))])

(with-expresso [+ * - /]
  (transform-with-rules group-rules (+ 1 0))  ;=> 1
  (transform-with-rules group-rules (* 3 1)) ;=> 3
  ;works because + and * are commutative
  (transform-with-rules group-rules (+ 0 1)) ;=> 1
  (transform-with-rules group-rules (* 1 3)) ;=> 3
  (transform-with-rules group-rules (+ 3 (- 3))) ;=> 0
  (transform-with-rules group-rules (* 3 (/ 1 3))) ;=> 1
  (transform-with-rules group-rules 
    (+ (* 1 (* 3 (/ 1 3))) (+ (* 1 2) (- (* 2 1))))) ;=> 1
)

3 What's next

In the next week I want to work further on the transformation with rules, try different strategies (postwalk vs prewalk) also in a core.logic context, and give proper support for n-ary expressions. I plan to do this by allowing a &* syntax in the rules which matches the rest of the arguments.

Ang again, please comment.

Date: 2013-06-23T13:36+0200

Author: Maik Schünemann

Org version 7.8.09 with Emacs version 23

Validate XHTML 1.0