Symbolic Matrix Library in clojure, advice needed

Question

While there are plenty of java libraries for linear algebra, Clojure currently does not have an idiomatic computer algebra system that would include support for symbolic math. As a start I figure I can start working on something simple.

As a first step, I think getting the data structures right would be a good start.

Step 1: Implement a persistent matrix

I will be using deftype (or reify), and for now, for ease of implementation, I will use a hashmap for storage (Please suggest an alternative if think its better but state tradeoffs). (One can imagine many different implementations depending on your performance requirements, such as using arrays or delegating to an external java library and implementing some sort of transients interface.)

My question is, what interfaces/protocols should I be considering to implement? (In general, what is a good listing of all the protocols/interfaces that clojure uses?) Also is there any advice on how to implement these?

My list of things to implement:

-Assoc'ing would be useful, to modify sections of the matrix in an immutable manner

-treating the matrix as a function as an accessor of the elements, I was thinking you could pass a two-tuple to return a single element, a single value (index by width*y+x), hashmap to get columns, rows, or minor, via custom query hashmap/language.

Note, my goal at the moment is to design good abstractions that enable flexibility in choosing implementations.

Answer 1

I more or less manage the linear algebra modules in SymPy, a major Python symbolics package. I'll give you my perspective coming from a traditional computer algebra system.

We have three separate implementations for three important use cases

1. Mutable matrices -- Despite being in Python SymPy is immutable by default. We've actually broken this rule for matrices. Matrix algorithms are the standard example of where you really need to switch to mutability for performance reasons.

2. Immutable matrices -- But you'd like to have the option to switch back. Our intended workflow is as follows

   1. Build an immutable matrix
   2. Switch to mutability and perform some algorithm
   3. Switch back to immutability and present this to the user

3. Matrix Symbols -- Often you don't need to deal with explicit entries in the matrix but would rather deal with the idea of a matrix. See this scicomp.stackexchange post. This is my current work and I find it to be very exciting.

There are other splits like dense representations vs sparsity. Symbolic linear algebra is a big and important field. I look forward to seeing the Clojure community's collective solution.

Answer 2

Although it is a Java library, I designed vectorz to be used from Clojure.

It offers a lot of data structures and algorithms for high performance vector and matrix maths. You might find it useful. I'm currently using it for both computer graphics and machine learning in Clojure.

The matrices and vectors are mutable, but I found that this was a necessary evil: using immutable vectors and matrices was simply too slow for many algorithms.

I'd be interested in building an idiomatic clojure wrapper (including immutable versions of vectors and matrices) if enough people would find this useful and/or want to get involved.

Answer 3

It may be of interest to readers of this thread that Gerry Sussman's scmutils system is being ported to Clojure. This is a very advanced CAS, offering things like automatic differentiation, literal functions, etc, much in the style of Maple. It is used at MIT for advanced programs on dynamics and differential geometry, and a fair bit of electrical engineering stuff. It is also the system used in Sussman&Wisdom's "sequel" (LOL) to SICP, SICM (Structure and Interpretation of Classical Mechanics). Although originally a Scheme program, this is not a direct translation, but a ground-up rewrite to take advantage of the best features of Clojure. It's been named sicmutils, both in honour of the original and of the book This superb effort is the work of Colin Smith and you can find it at https://github.com/littleredcomputer/sicmutils .

I believe that this could form the basis of an amazing Computer Algebra System for Clojure, competitive with anything else available. Although it is quite a huge beast, as you can imagine, and tons of stuff remains to be ported, the basics are pretty much there, the system will differentiate, and handle literals and literal functions pretty well. It is a work in progress. The system also uses the "generic" approach advocated by Sussman, whereby operations can be applied to functions, creating a great abstraction that simplifies notation no end.

Here's a taster:

> (def unity (+ (square sin) (square cos)))
> (unity 2.0)  ==>  1.0
> (unity 'x)   ==> 1 ;; yes we can deal with symbols
> (def zero (D unity))  ;; Let's differentiate
> (zero 2.0)   ==> 0

SicmUtils introduces two new vector types “up” and “down” (called “structures”), they work pretty much as you would expect vectors to, but have some special mathematical (covariant, contravariant) properties, and also some programming properties, in that they are executable!

> (def fnvec (up sin cos tan))  => fnvec
> (fnvec 1)   ==> (up 0.8414709848078965 0.5403023058681398 1.5574077246549023)
> ;; differentiated
> ((D fnvec) 1)  ==>  (up 0.5403023058681398 -0.8414709848078965 3.425518820814759) 
> ;; derivative with symbolic argument
> ((D fnvec) 'θ) ==> (up (cos θ) (* -1 (sin θ)) (/ 1 (expt (cos θ) 2)))  

Partial differentiation is fully supported

> (defn ff [x y] (* (expt x 3)(expt y 5)))
> ((D ff) 'x 'y) ==> (down (* 3 (expt x 2) (expt y 5)) (* 5 (expt x 3) (expt y 4))) 
> ;; i.e. vector of results wrt to both variables

The system also supports TeX output, polynomial factorization, and a host of other goodies. Lots of stuff, however, that could be easily implemented has not been done purely out of lack of human resources. Graphic output and a "notepad/worksheet" interface (using Clojure's Gorilla) are also being worked on.

I hope this has gone some way towards whetting your appetite enough to visit the site and give it a whirl. You don't even need Clojure, you could run it off the provided jar file.

=========

PS. Incidentally, to answer the original questiin directly, yes, sicmutils does support symbolic structures: you can set up a matrix representation where the entries are formulae, e.g. a matrix of rotation, and then evaluate (multiply) it for a given coordinate. It's wonderfully flexible in that way.

Answer 4

This may not be answering your question, but one thing I found while working with Incanter was I needed to be able to access elements with wraparound. Same for vectors, it's sometimes handy to be able to provide negative indices to access elements in reverse from the end, or out-of-range indices to access elements from the start. Naturally computing the offsets adds overhead but it's sometimes a feature you're willing to pay for

Answer 5

MPL is a simple symbolic mathematics library written in portable R6RS Scheme. Here's a short introduction.

Since Scheme is a Lisp, MPL should be quite straightforward to port to Clojure.