Performance tuning for a simple SLE solver
andi.x...@googlemail.com
I'm new to Clojure, and let me first say that I love it! At least I
love the language, but I have some concerns regarding performance:
My first try was to implement a Gauß elemination algorithm for solving
a system of linear equations. Here is the code:
http://www.xenoage.com/extern/zongblog/matrix.clj
Solving 1.000.000 SLEs (see last line) took 33.000 ms on my machine,
while the corresponding algorithm written in Java needed less than 300
ms.
Since I am a Java programmer and I have no experience in functional
programming, I probably have to apologize for the bad code. Anyway, I
already tried to make it faster using the well-known performance tips
(like type hints), but they did not really help (ok, 25.000 ms instead
of 33.000, but it is still too much).
What are my options?
- Can you identify problems in my code? I do not mean things that make
it 10% faster, but at least double as fast.
- Could you say, that Clojure is not made for such numerical things,
and I should use Java code for these performance-critical algorithms?
(would be perfectly ok for me)
- Should I try it in Clojure, but using a (mutable) Java array instead
(also ok for me, since this mutable structure is only used locally and
within a single thread)
Any ideas? Thanks in advance :-)
Andi
B Smith-Mannschott
(I'm still a Clojure noob too, so...)
This is the kind of code, I would just write in Java.
However: your Clojure implementation makes use of exact arithmetic,
i.e. rational numbers. Does the Java version do that too, or is it
using floating point? If so, then the two benchmarks are not directly
comparable.
I hacked up your implementation to use http://clojure.org/transients
which were added to Clojure after the 1.0 release. This helped
performance some. See comments in code below.
// Ben
--------------------------------------
(ns matrixt)
; This module applies the Gauss algorithm to a matrix to solve a SLE.
;
; @author Andreas Wenger
; (hacked up to use transients -- bsmit...@gmail.com)
;; Used transients instead of normal persistent datastructures within
;; solver since transients are more efficient to update in place at
;; the cost of invalidating previous versions of themselves. (thus not
;; persistent.)
;;
;; This has helped performance, but nothing earth-shattering. Haven't
;; tried type hinting -- or indeed profiling. Also, this is the first
;; time I've tried to use transients for something and I find the results
;; rather ugly. (tvec-map! especially)
;;
;; On a 1.6 GHz Atom (Dell Mini 9)
;; (do (time (dotimes [_ 100000] (matrixt/solve matrixt/testsle2)))
;; (time (dotimes [_ 100000] (matrix/solve matrix/testsle2))))
;; "Elapsed time: 12816.971679 msecs"
;; "Elapsed time: 18283.428481 msecs"
; matrix to solve:
; AAA b
; AAA b
; AAA b
(defstruct sle :A :b)
(def testsle1 (struct sle [[3 2][2 4]] [1 7]))
(def testsle2 (struct sle [[1 2 3][2 5 7][8 2 1]] [1 7 3])) ;
solution: [-16/9 110/9 -65/9]
(def testsle3 (struct sle [[1 2 3][0 1 3][0 0 1]] [2 4 6]))
(defn tvec-map!
"Map function over elements of tvec (a transient vector). Replace
elements from tvec with corresponding result computed by f."
([f tvec]
(let [n (count tvec)]
(loop [i 0, tvec tvec]
(if (< i n)
(recur (inc i) (assoc! tvec i (f (nth tvec i))))
tvec))))
([f tvec1 tvec2]
(let [n (count tvec1)]
(loop [i 0, tvec tvec1]
(if (< i n)
(recur (inc i) (assoc! tvec i (f (nth tvec i) (nth tvec2 i))))
tvec)))))
(defn sle-transient [s]
"Transform SLE structmap into a transient map containing transient
vectors for faster inplace updating."
(transient {:A (transient (vec (map transient (:A s))))
:b (transient (:b s))}))
(defn sle-persistent! [s]
(struct sle
(persistent! (tvec-map! persistent! (:A s)))
(persistent! (:b s))))
(defn normalize-line!
"Converts the line of a given SLE (parameter s) with the given
index (parameter l) into the form that the column with the given
index has a 1."
[s l]
(let [dividend (get-in s [:A l l])
div #(/ % dividend)]
(-> s
(assoc! :A (assoc! (:A s) l (tvec-map! div (get-in s [:A l]))))
(assoc! :b (assoc! (:b s) l (div (get-in s [:b l])))))))
(defn subtract-line!
"Subtracts the line with the given index (ldest) of a given
SLE (parameter s) by an appropriate multiple of the line with index
lsrc"
[s lsrc ldest]
(let [factor (get-in s [:A ldest lsrc])
minus #(- %1 (* factor %2))]
(-> s
(assoc! :A (assoc! (:A s) ldest (tvec-map! minus
(get-in s [:A ldest])
(get-in s [:A lsrc]))))
(assoc! :b (assoc! (:b s) ldest (minus (get-in s [:b ldest])
(get-in s [:b lsrc])))))))
(defn forward-lines!
"Converts the lines of a given SLE (parameter s) after the given
index (parameter l) into a form with the first l columns containing
0. This is done by subtracting an appropriate multiple of the line
with index l"
[s l]
(let [n (count (:A s))]
(loop [acc s, i (inc l)]
(if (>= i n)
acc
(recur (subtract-line! acc l i) (inc i))))))
(defn forward-lines
[s l]
(sle-persistent! (forward-lines! (sle-transient s) l)))
(defn forward!
"Converts the given SLE (parameter s) into a triangular form"
[s]
(let [n (count (:A s))]
(loop [acc s, i 0]
(if (>= i n)
(normalize-line! acc (dec n))
(recur (forward-lines! (normalize-line! acc i) i) (inc i))))))
(defn backward-lines!
"Converts the lines of a given SLE (parameter s) before the given
index (parameter l) into a form with the last n-l columns
containing 0 (n is the number of variables in the SLE). This is
done by subtracting an appropriate multiple of the line with index l"
[s l]
(loop [acc s, i (dec l)]
(if (< i 0)
acc
(recur (subtract-line! acc l i) (dec i)))))
(defn backward!
"Backward substitution"
[s]
(let [n (count (:A s))]
(loop [acc s, i (dec n)]
(if (< i 1)
acc
(recur (backward-lines! acc i) (dec i))))))
(defn solve
"Solves the given SLE"
[s]
(-> s sle-transient forward! backward! sle-persistent! :b))
;(time (dotimes [_ 1000000] (solve testsle2)))
luskwater
them to nested calls on the map or vector:
(my-get-in sle [:A 1 2]) =>
(((sle :A) 1) 2)
The form (data key) seems to be about 10 to 15 times faster than the
equivalent (get-in data [key]) on my laptop.
It cut the runtime in half or so. I would think that floating point
operations might be faster, too.
On Oct 9, 6:38 am, "andi.xeno...@googlemail.com"