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 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))
  (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)
(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

Freitag, 14. Juni 2013

Constructing Algebraic Expressions

Constructing Expressions

Constructing Expressions

Gsoc finally starts monday. And for the start I want to have a good idea of how the user interface of expresso should look like, for this I want to discuss the best way of aktually constructing expresso expressions :).

1 Representation of expresso expressions

Because core.logic currently doesn't play nicely with custom types and algebraic expressions fit nicely with normal clojure s-expressions I plan to represent algebraic expressions just as that. However, that doesn't mean it is straightforward to construct them - expresso has to know the full qualified function symbol to differentiate between - for example - clojure.core/* and the core.matrix * which behave differently. Also other mathematical properties should be stored about the functions - I plan to use metadata for it.

2 How should using expresso feel like

Hopefully as seamless as possible.

One approach would be a macro - create-expression. It would construct the enhanced form of the s-expression from the s-expression it gets. This, however has some drawbacks. First the user either needs to use ~ to aktually get the data in the expressions. Or you have to have map as argument which would map from symbols to values. Further than that, because it is a macro, it can't be used in higher order functions.

(let [matrix [[1 0 1]
              [0 1 0]
              [1 0 1]]]
  (simplify (construct-expression (* 5 (determinant ~matrix))))
  (simplify (construct-expression (* 5 (determinant x)) {'x matrix})))

The other approach would be to have construction functions in an own namespace, with concise names, just like core.matrix did with the operators namespace.

(ns test
  (:use [])
  (:refer-clojure :exclude [* - + == /]))
(let [matrix [[1 0 1]
              [0 1 0]
              [1 0 1]]]
 (simplify (* 5 (determinant matrix))))

Mike Anderson also posted some interesting usage examples (with a possible syntax) on the clojure mailing list:

; 1. Solving simultaneous equations expressed in arrays:

(solve [x y] (== [x y] (* [[3 0] [0 10]] [y 1])))
=> [30 10]

;2a. Derivatives of array expressions:

(derivative [x] (* [x 2] [3 (* x x)]))
=> [3 (* 4 x)]

;2b. Derivatives that are themselves vectors / arrays:

(derivative [[x y]] (* [x 2] [3 (* x x)]))
=> [[3 (* 4 x)] 0]

;3. Simplification of computations on arrays:

(simplify (trace [[x y] [z (- 5 x)]]))
=> 5

;4. Simplification of array-level expressions

(simplify (* 10 A (some-identity-matrix) (inverse A) B))
=> (* 10 B)

What would be your favorite way of using expresso?

3 expresso on core.matrix expressions

Considering the close relation between expresso and core.matrix it is very likely that expresso (especially the optimizer) will be used on core.matrix expressions. Mike's examples also point in this direction. What about making expresso an implementation of core.matrix ? It would not calculate the results of matrix operations, but create the expressions under the hood. That way one could just switch the implementation of core.matrix to make the matrix expression available to expresso.

For the core.matrix folks: Do you think that is possible? There are certainly a few caveats in it but it seems to be the perfect way of using expresso in a core.matrix context.

Date: 2013-06-14T23:09+0200

Author: Maik Schünemann

Org version 7.8.09 with Emacs version 23

Validate XHTML 1.0

Mittwoch, 29. Mai 2013

GSOC Project - Algebraic Expressions

GSOC Project - Algebraic Expressions

GSOC Project - Algebraic Expressions

Hello clojurians. I am glad to say that my Proposal for the Algebraic Expression project idea was accepted! It is a very interesting project and can be a big benefit for the scientific computing interested clojure programmers. In this blog post I want to tell what I plan to achieve in expresso ( a good name for the project :)), and also to discuss some relevant design decisions I have to make. I really appreciate every comment or advice. The more I get, the better I can make the library for you! So please comment.

1 General Idea

Expresso aims to be a general purpose library to syntactically manipulate algebraic expressions, and combines clojure and the declarative power of core.logic to do this.

It will contain a solver and an optimizer, which are inspired by sympy and theano.

To give a more concrete idea, here ist the synopsis part of my proposal:

The goal is to implement a library for manipulating and analyzing algebraic expressions. This would allow several advanced features like solving a system of equations in regard to a variable X, or optimizing an expression for executing via core.matrix. The library should also be extensible, to build a foundation for symbolic manipulation of expressions, on which other more sophisticated transformations can be based.

I plan to use the search and backtracking capabilities of core.logic - the clojure port of mini-kanren - to manipulate expressions. The scope of the library would be:

  • setting up a logic programming context, in which manipulations can be applied
  • providing a general rule-based transformation facility - with convenient syntax for defining rules. A rule can be any logical clause, which enables even whole manipulations to be stated as rules
  • transforming an expression into a (or more) suitable standart form (i.e unnesting sums, eliminating zero multiplication,…) for further manipulations.
  • an input namespace, which makes it easy to transform and existing expression using core.matrix to a suitable one for algebraic manipulation

Also, it should contain two more sophisticated transformations

  • an optimizer, which transforms an expression to be used with core.matrix to an equivalent expression, which runs faster. It could do optimizations, like removing redundant calculations,moving computation to compile time, analyzing the expression for parrallel executiion, …
  • a solver for solving (sets of) equations. It would analyze the expression and recursively try to e.g isolate the variable if it occurs only once, solve the equation as an polynominal, contract multiple occurences of the variable. each of this transformations above can be a rule by itself and can contain rules by itself. This would make the solver extensible.

2 Core design decisions

Here are the core design decisions I need to address in order to be able to implement the more sophisticated solver and optimizer on top of it

  • how to encode an expression
    I currently plan to encode expressions as just plain old s-expressions. An extra type will be problematic because s-expressions are easy to handle with core.logic. I did not manage to get a self defined type play well in all situations with core.logic, but maybe some core.logic specialist could give advice if it is possible and recommended to use a custom type with core.logic. This raises other issues such as
    • core/* and core.matrix/* behave differently - symbols will have to be fully qualified but there are many types on which + behaves the same - So I also want to have one set of rules for these +'ses. Clojures ad-hoc hierarchies could help here.
    • I will need to store additional information on the expressions, like mathematical properties wich can be used to optimize it. without a custom type I could use metadata for it.
  • how to express a rewrite rule
    expresso will heavily rely on rewrite rules and also let the set of rules be extensible by users
    • I plan to use a core.logic clause as a rule. It will succeed when the rule is applicable and fail when not when succeeding it unificates the second argument with the simplified rule.
    • this would allow arbitrary complex rules, such as a whole calculation routing as a rule
    • It will be useful to have a good syntax to create rules - I currently use a very simple one based on matche
  • how to apply rules on an expression
    • there are many ways to traverse an expression tree, it is probably best to support multiple and see which performs best
    • rule based translators tend to perform badly with many rules. There need to be a good strategie to reduce the time to search for the right rule.
  • how existing rule-based translators are extensible by users
    • where to add a new rewrite rule - order matters
  • how to construct an expression
    Use of the library should be so seamless as possible. It would be perfect, if it could mirror the core.matrix api in an input namespace to construct an expression on which the library can work.

This are the central issues I can think of for now. There are others for the solver and optimizer, and surely other for the core part which I forgot…

3 The code

Mike Anderson created the expresso library and I committed to it to test my ideas. For now it contains a proof of concept rule-based translator and solver. In the code snippet below there are the rules used by simplification and solving and test cases. Remember that it is only a proof of concept.

  ; a rule to calculate an expression - is just a simple core.logic clause
(defn calculo [eq eq1]
  (no-variablo eq)
  (resulto eq eq1))

; A rule is a clause which gets an expression and returns a modified expression
; If a rule is not applicable, it fails
(def simp-rules
  [(rule ['+ 0 x] :=> x)
   (rule ['+ x 0] :=> x)
   (rule ['* 0 x] :=> 0)
   (rule ['* x 0] :=> 0)
   (rule ['* 1 x] :=> x)
   (rule ['* x 1] :=> x)
   (rule ['- 0 x] :=> x)
   (rule ['- x 0] :=> x)
   (rule ['- x x] :=> 0)
   (rule ['* y ['+ 'x z]] :=> ['+ ['* y 'x] ['* y z]])
   (rule ['* ['+ 'x z] y] :=> ['+ ['* y 'x] ['* y z]])
   (rule ['+ x x] :=> ['* 2 x])
   (rule ['+ ['* y 'x] ['* z 'x]] :=> ['* ['+ z y] 'x])
   (rule ['+ ['* y 'x] ['* 'x z]] :=> ['* ['+ z y] 'x])
   (rule ['+ ['* 'x y] ['* z 'x]] :=> ['* ['+ z y] 'x])
   (rule ['+ ['* 'x y] ['* 'x z]] :=> ['* ['+ z y] 'x])

 ;simple example
  (run* [q] (simplifyo '(+ x 0) q))
  ;=> (x)
  ;rules are applied recursively
  (run* [q] (simplifyo '(* x (+ 0 (* 3 (* x 0)))) q))
  ;=> (0)

  ;calculate is also just a rule, which checks if there is no variable in expression and calculates it
  (run* [q] (simplifyo '(+ (* 3 4) (+ (* x 0) 7)) q))
  ;=> (19)
; will be used from the solver 1 means the variable x is on the left side
(def isolate-rules
   (rule [1 ['= ['+ x y] z]] :=> ['= x ['- z y]])
   (rule [2 ['= ['+ y x] z]] :=> ['= x ['- z y]])
   (rule [1 ['= ['* x y] z]] :=> ['= x ['/ z y]])
   (rule [2 ['= ['* y x] z]] :=> ['= x ['/ z y]])
   (rule [1 ['= ['- x y] z]] :=> ['= x ['+ z y]])
   (rule [2 ['= ['- y x] z]] :=> ['= ['- 0 x] ['- z y]])

(run* [q] (solveo '(= (+ (* 3 (+ (* x 4) (* x 3))) 5) 89) 'x  q))
 ;=> ('(= x 4))

You can see the code here. The wiki pages contain more information. If you have any feature request or comment/advice for me, please comment or post on the mailing list. If you already have a use case on which you would want to use the library please tell me so I can have an accurate picture on what forms of expressions it should be able to handle and how it can be best used.

I am really glad that I managed to get into gsoc with such a good project. I am sure the summer will be interesting!

Date: 2013-05-29T16:53+0200

Author: Maik Schünemann

Org version 7.8.09 with Emacs version 23

Validate XHTML 1.0