Packing and Unpacking Numbers : Deterministic Randomness

;; Packing and Unpacking Numbers

;; This post is kind of pointless, since while writing it I discovered a better way to do the same thing.

;; But I also think it's kind of cool, so I'm putting it up in case anyone else finds a use for it.

;; Suppose we've got a function which gives random results (I mean, deliberately...)

(defn randomfunction []
  (list 
   (case (rand-int 2) 0 '+ 1 '-)
   (rand-int 100)
   (list 
    (case (rand-int 2) 0 '+ 1 '-)
    (rand-int 100)
    (rand-int 100))))

(randomfunction) ;-> (+ 43 (- 71 75))
(randomfunction) ;-> (- 18 (- 13 22))
(randomfunction) ;-> (- 36 (+ 60 3))
(randomfunction) ;-> (+ 56 (+ 34 4))
(randomfunction) ;-> (- 87 (- 61 76))

;; And it comes to testing it:

(= (randomfunction) '(+ 26 (- 69 0))) ;-> false

;; We'd prefer it to be deterministic. 

;; So we can separate the concerns thus:

(defn deterministicfunction [a b c d e]
  (list 
   (case a 0 '+ 1 '-)
   b
   (list 
    (case c 0 '+ 1 '-)
    d
    e)))

(defn randomfunction []
  (deterministicfunction (rand-int 2) (rand-int 100) (rand-int 2) (rand-int 100) (rand-int 100)))

;; tadaa!

(= (deterministicfunction 0 26 1 69 0) '(+ 26 (- 69 0))) ;-> true

(randomfunction) ;-> (+ 30 (+ 52 95))
(randomfunction) ;-> (- 97 (- 57 61))
(randomfunction) ;-> (+ 88 (+ 26 15))

;; Although useful, this is a right royal pain in the neck, and horribly error-prone

;; It occurred to me that it might be made less painful by setting the
;; seed of the random number generator before running tests, but I
;; also had a need to make 'deterministic random' values in the code
;; itself, and having everything hanging off a global variable that
;; any old function occasionally resets struck me as a recipe for
;; disaster.

;; So I thought, well why not pack all the random parameters into one big number:

(def << bit-shift-left)

;; If all the random ranges are powers of two then this is easy:

(let [[a b c d] [(rand-int 8) (rand-int 4) (rand-int 8) (rand-int 2)]]
 [[a b c d] (+ (<< a (+ 1 3 2)) (<< b (+ 1 3)) (<< c 1) (<< d 0))]) 
; [[2 2 4 0] 168] 
; [[3 1 4 0] 216] 
; [[5 3 6 0] 380]

;; Looking at the last one
;; 5/3/6/0 is 
;; 101/11/110/0 in binary

;; and binary 101111100 is:
2r101111100 ; 380

;; in decimal

;; To extract the numbers, we can do:

(mod 380 2) ; 0 
(/ (- 380 (mod 380 2)) 2) ; 190 
(mod 190 8) ; 6 
(/ (- 190 (mod 190 8)) 8) ; 23
(mod 23 4) ; 3
(/ (- 23 (mod 23 4)) 4) ; 5
(mod 5 8) ; 5

;; And this motivates the following recursion:

(defn ^:dynamic unpack [rl n]
  (if (empty? rl) '() 
      (let [[rh & rr] rl]
        (cons (mod n rh) (unpack rr (/ (- n (mod n rh)) rh))))))

;; Which behaves as I was hoping

(unpack '() 0) ;-> ()
(unpack '() 1) ;-> ()
(unpack '(8) 1) ;-> (1)
(unpack '(8) 9) ;-> (1)
(unpack '(8 8) 9) ;-> (1 1)

(unpack '(2 8 4 8) 380) ;-> (0 6 3 5)
(unpack '(2 8 4 8) 216) ;-> (0 4 1 3)
(unpack '(2 8 4 8) 168) ;-> (0 4 2 2)

;; And that then motivates the following packing function

(defn ^:dynamic pack [rlist vlist]
  (if (empty? rlist) 0
      (let [[rh & rr] rlist
            [vh & vr] vlist] 
        (+ (* rh (pack rr vr)) vh))))

;; Which also works
(pack '()()) ;-> 0
(pack '(8) '(1)) ;-> 1
(pack '(8 8) '(1 1)) ;-> 9
(pack '(2 8 4 8) '(0 6 3 5)) ;-> 380
(pack '(2 8 4 8) '(0 4 1 3)) ;-> 216
(pack '(2 8 4 8) '(0 4 2 2)) ;-> 168

;; But the cool thing is that it works whether or not the ranges are powers of two:

(use 'clojure.tools.trace)

(let [a (for [i (range (inc (rand-int 7)))] (rand-int 30))
      b (for [i a] (rand-int i))]
  (clojure.tools.trace/dotrace [unpack pack] (unpack a (pack a b))))

;; TRACE t3191: (pack (12 23 5 5 8 3) (4 9 3 4 0 0))
;; TRACE t3192: | (pack (23 5 5 8 3) (9 3 4 0 0))
;; TRACE t3193: | | (pack (5 5 8 3) (3 4 0 0))
;; TRACE t3194: | | | (pack (5 8 3) (4 0 0))
;; TRACE t3195: | | | | (pack (8 3) (0 0))
;; TRACE t3196: | | | | | (pack (3) (0))
;; TRACE t3197: | | | | | | (pack nil nil)
;; TRACE t3197: | | | | | | => 0
;; TRACE t3196: | | | | | => 0
;; TRACE t3195: | | | | => 0
;; TRACE t3194: | | | => 4
;; TRACE t3193: | | => 23
;; TRACE t3192: | => 538
;; TRACE t3191: => 6460
;; TRACE t3198: (unpack (12 23 5 5 8 3) 6460)
;; TRACE t3199: | (unpack (23 5 5 8 3) 538)
;; TRACE t3200: | | (unpack (5 5 8 3) 23)
;; TRACE t3201: | | | (unpack (5 8 3) 4)
;; TRACE t3202: | | | | (unpack (8 3) 0)
;; TRACE t3203: | | | | | (unpack (3) 0)
;; TRACE t3204: | | | | | | (unpack nil 0)
;; TRACE t3204: | | | | | | => ()
;; TRACE t3203: | | | | | => (0)
;; TRACE t3202: | | | | => (0 0)
;; TRACE t3201: | | | => (4 0 0)
;; TRACE t3200: | | => (3 4 0 0)
;; TRACE t3199: | => (9 3 4 0 0)
;; TRACE t3198: => (4 9 3 4 0 0)
;; (4 9 3 4 0 0)

;; Which I found slightly smug-making

;; So now we can say:
(defn deterministicfunction [big-rand]
  (let [[a b c d e] (unpack [2 100 1 100 100] big-rand)]
    (list 
     (case a 0 '+ 1 '-)
     b
     (list 
      (case c 0 '+ 1 '-)
      d
      e))))

(defn randomfunction []
  (deterministicfunction (rand-int (* 2 100 2 100 100))))

;; Behold the many benefits of this scheme:

;; We can enumerate all possible varieties
(deterministicfunction 1) ;-> (- 0 (+ 0 0))
(deterministicfunction 2) ;-> (+ 1 (+ 0 0))
(deterministicfunction 3) ;-> (- 1 (+ 0 0))

;; We can generate examples using a single number
(deterministicfunction 10345778979) ;-> (- 89 (+ 94 88))

;; Or we can get random examples
(randomfunction) ;-> (+ 99 (+ 72 8))

;; If we are not desperately fussed about the statistical properties of things
;; we do not need to worry about the exact range of the random number at all:
;; I think that this will produce only tiny biases:

(defn randomfunction []
  (deterministicfunction (rand-int (<< 1 30))))

;; And it still looks pretty random to me
(randomfunction) ;-> (+ 76 (+ 10 87))
(randomfunction) ;-> (+ 76 (+ 81 3))
(randomfunction) ;-> (- 27 (+ 18 63))

Clojure Emacs Eval Paste Keyboard Macros / Generating your Regression Tests Automatically while Writing your Functions

;; Clojure Emacs Eval Paste Keyboard Macros / 
;; Generating your Regression Tests Automatically while Writing your Functions

;; There's a phenomenally useful feature of emacs when using clojure
;; and nrepl where you can evaluate the result of an
;; expression and paste the result directly into a buffer.

;; However my attempts to use it in keyboard macros have always led to frustration.

;; I asked about it on stack overflow:
;; http://stackoverflow.com/questions/14959155/how-can-i-make-emacs-keyboard-macros-work-properly-when-pasting-the-results-of-c#14960321

;; and a nice chap called sds helped me come up with the following work-around:

;; Say I've defined the function:

(defn << [a b]
  (if (zero? b) a
      (<< (* 2 a) (dec b))))


;; And now I want to check that it works:

;; I evaluate the following two forms in emacs lisp (in the *scratch* buffer)

(defun clojure-eval-paste-test ()
  (interactive)
  (next-line)
  (move-beginning-of-line nil)
  (insert "(is= ")
  (forward-sexp)
  (insert " ")
  (nrepl-eval-last-expression 't)
  (sleep-for 0.1)
  (insert ")"))

(global-set-key [f5] #'clojure-eval-paste-test)


;; Now I type:

(<< 1 0)
(<< 1 1)
(<< 1 2)
(<< 2 0)
(<< 2 2)
(<< 2 3)

;; And then I move point to just above the first line, and press f5 six times, and I get:


(is= (<< 1 0) 1)
(is= (<< 1 1) 2)
(is= (<< 1 2) 4)
(is= (<< 2 0) 2)
(is= (<< 2 2) 8)
(is= (<< 2 3) 16)

;; Which both reassures me that it does work, and provides the skeleton of a regression test for it!



;; There are more elegant ways of doing this, presumably, one being to
;; do as sds suggests, and code up the function using lower level
;; calls rather than the (sleep-for 0.1) thing that will obviously go
;; wrong if the eval takes too long.

;; And another way would be to fix it so that nrepl-eval-last-expression actually works in keyboard macros. 

;; But this will do for now, and gets rid of a papercut bug that's been annoying me for ages.

Mathematics as an Experimental Science

;; Mathematics as an experimental science

;; Beta-Bernoulli distribution / Machine Learning

;; I'm reading Kevin Murphy's Machine Learning: A Probabilistic Perspective.

;; As one of the exercises in chapter three, I've just proved that:

;; "For the Beta-Bernoulli model, the probability of the data is the
;; ratio of the normalizing constants of the prior distribution and
;; the posterior distribution". 

;; What I think this actually means is:

;; If you start with Coins which are biased, and the bias is chosen
;; uniformly (so that for instance a fair coin is just as likely as a
;; coin which comes up heads 1/3 of the time, or a coin which always
;; comes up heads).

;; And you pick out the fair ones by tossing them all 20 times, and
;; throw away all the ones that don't come up heads exactly 10 times
;; and tails 10 times.

;; Then, because your test has not quite guaranteed you fair coins,
;; but only coins which are quite a bit fairer than you started with:

;; You are ever so slightly more likely to see three heads (or three
;; tails) in three tosses than you would be if the coins were truly
;; fair.

;; In fact the bias is now the same as if the coins had been chosen
;; from the Beta(11,11) distribution.

;; And so the chances of getting three heads in a row is not, as you
;; might have naively expected, one in eight,

;; But rather B(14,11)/B(11,11), where B is the Beta function

;; And that works out to be G(14)G(11)G(22)/G(25)G(11)G(11), where G is the Gamma function.

;; And that is 13!10!21!/24!10!10!, where ! is the factorial function

;; And that is (11*12*13)/(22*23*24), where / and * need no introduction if you have got this far in the post.

;; And that is:

(/ (* 11 12 13) 22 23 24) ;-> 13/92

;; or 

(float (/ 13 92)) ; 0.14130434

;; Which is to say, slightly more than one in eight.

;; And I am inordinately pleased with myself, not for having proved
;; this result, which was easy, but for having worked out what the
;; mysterious squiggles in the book might actually mean in practice.

;; And it occurs to me that Mathematics is in fact an experimental
;; science, which makes definite predictions about physical things,
;; such as the movements of electrons in my computer when I type the
;; following things:

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;; Let us make coins with biases distributed uniformly between 0 and 1:

(defn make-coin []
  (let [theta (rand)]
    (fn [] (if (< (rand) theta) :h :t))))

;; Here is one:
(def a (make-coin))

;; Let us test it
(frequencies (for [i (range 20)] (a))) ;-> {:h 11, :t 9}

;; It is no good. To the scrapheap with it:
(ns-unmap *ns* 'a)

;; Let us instead make a large number of coins:
(def coins (for [i (range)] (make-coin)))

;; And throw away all the ones that do not satisfy our criterion
(def fair-coins
  (filter (fn[coin] (= {:t 10 :h 10 } (frequencies (for [i (range 20)] (coin))))) coins))

;; Let us toss each of our fair coins three times:
(def results (for [coin fair-coins] (frequencies (for [i (range 3)] (coin)))))

;; And add up all the results
(def collected-results (drop 1 (reductions (fn[m f] (assoc m f (inc (get m f 0)))) {} results)))

;; And calculate the empirical distribution of those results:
(def empirical-distribution (map 
                             (fn[m] [(float (/ (m {:h 3} 0) (reduce + (vals m))))
                                     (float (/ (m {:t 3} 0) (reduce + (vals m))))])
                             collected-results))

;; It does seem to me that after a while, these numbers settle down to something near 14%
(doseq [e empirical-distribution] (println e))

A Naive Bayesian Classifier

Anyway this very evening Hermann Hauser of this parish is telling me that Machine Learning is the coming next big thing and those skills are highly in demand.

So if anybody wants to give me a fuckload of money for this old rope, they must feel free

;; A Naive Bayesian Classifier

;; Ed Jackson ( http://boss-level.com ) and I are currently working
;; our way through Kevin Murphy's book:
;; Machine Learning: A Probabilistic Perspective.

;; Ed tells me that there is big money in this Machine Learning game.

;; We have just read chapter 3, which involves a hair-raising amount
;; of Greek, but when all is said and done boils down to this:

;; We got Martians, who tend to be thin, tall and green
(defn make-martian []
  #{:martian
    (if (< (rand) 0.2) :fat :thin)
    (if (< (rand) 0.7) :tall :short)
    (if (< (rand) 0.8) :green :blue)})

; And we got Venusians, who tend to be fat and blue
(defn make-venusian []
  #{:venusian
    (if (< (rand) 0.8) :fat :thin)
    (if (< (rand) 0.5) :tall :short)
    (if (< (rand) 0.3) :green :blue)})

; They visit Harry's 
(defn make-space-bar-patron []
  (if (< (rand) 0.3) (make-martian) (make-venusian)))

; Harry takes notes
(def client-database
  (for [i (range 10000)] (make-space-bar-patron)))

; Like this:
(take 5 client-database)
(#{:tall :martian :blue :thin} 
 #{:blue :venusian :short :fat} 
 #{:blue :venusian :short :thin} 
 #{:tall :martian :green :fat} 
 #{:martian :blue :short :thin})

; And tabulates them so:
(def martians (filter :martian client-database))
(def venusians (filter :venusian client-database))

(count martians) ; 3054
(count (filter :thin martians)) ; 2454
(count (filter :green martians)) ; 2438
(count (filter :tall martians)) ; 2177


(count venusians) ; 6946
(count (filter :fat venusians)) ; 5628
(count (filter :blue venusians)) ; 4919
(count (filter :tall venusians)) ; 3535


; A guy walks into Harry's
(def new-guy (make-space-bar-patron))

; But he is not wearing his uniform
(reduce disj new-guy [:martian :venusian])  ; #{:blue :short :fat}

; And so a sweepstake comes into being:
(defn probability [characteristic given-class prior-for prior-against]
  (let [class (if given-class (filter given-class client-database) client-database)
        for (count (filter characteristic class))
        against (- (count class) for)
        total-for (+ prior-for for)
        total-against (+ prior-against against)]
    (/ total-for (+ total-for total-against))))



(def m (* (probability :martian nil 1 1)      ; most people come in here are venusians
          (probability :blue  :martian 1 1)   ; and martians tend to be green
          (probability :short :martian 1 1)   ; also he's short. Martians are tall, often
          (probability :fat   :martian 1 1))) ; and fat. All the fat guys is from Venusburg

(def v (* 
        (probability :venusian nil 1 1)       ; most likely he's from venus
        (probability :blue  :venusian 1 1)    ; the blueness bears this out
        (probability :short :venusian 1 1)    ; sometimes they are short, sometimes they are tall
        (probability :fat   :venusian 1 1)))  ; but they are mostly fat

(float (/ m (+ m v))) ; 0.017495114
(float (/ v (+ m v))) ; 0.9825049

;; 45:1 he's from Venus !

new-guy ; #{:blue :short :venusian :fat}

;; As indeed he is.

Clojure is Fast / Is Clojure Still Fast? / Clojure is as Fast as Java

;; Is Clojure Still Fast ?

;; Once upon a time I wrote a blog post saying that clojure was fast.  It still is, and optimizing
;; it is now much easier than it used to be, but it doesn't seem to be *quite* as fast as it once
;; was.

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;; I optimized a program to solve a differential equation in a
;; simple-minded way:

;; The equation was: dy/dt = f(t,y) where f(t, y) = t-y and y=0 when t=0
;; ( the exact solution is: y=e^(-t)+t-1 )

;; Here's a program to solve it using Euler's method

(defn f [t y] (- t y))

(defn solveit [t0 y0 h its]
  (if (> its 0) 
    (let [t1 (+ t0 h)
          y1 (+ y0 (* h (f t0 y0)))]
      (recur t1 y1 h (dec its)))
    [t0 y0 h its]))

;; For every go round the loop we have to:
;; compare its with 0,
;; branch depending on the result,
;; add t0 to h,
;; call f with t0 and y0,
;; multiply h and the result,
;; add that to y0,
;; jump.

;; So if this was an assembly language program that worked the way
;; you'd expect, each loop would take 7 cycles.

;; According to /proc/cpuinfo, my netbook runs at 1.662 GHz.

(def *cpuspeed* 1.662)

;; We care about how many cycles each iteration of the solver takes:

(defmacro cyclesperit [expr its]
  `(let [start# (. System (nanoTime))
         ret# ( ~@expr (/ 1.0 ~its) ~its )
         finish# (. System (nanoTime))]
     (int (/ (* *cpuspeed* (- finish# start#)) ~its))))


;; With the program as written, this estimate turns out to have been a little optimistic.

;; The figures in the original post were on a desktop machine that I no longer have, which was more
;; powerful than my current netbook.

;; To two significant figures, the results of the timing expression
(cyclesperit (solveit 0.0 1.0) 1000000) 
;; are:

;; On my old desktop with clojure 1.2 : 2300 cycles

;; On my netbook with clojure 1.2 : 2800 cycles
;; On my netbook with clojure 1.3 : 2500 cycles
;; On my netbook with clojure 1.4 : 2400 cycles

;; So it looks like my netbook is not only slower in clock speed than my desktop was, but also in
;; terms of cycles/iteration. That's not surprising as the netbook has an Atom processor, optimized
;; for low power rather than for speed.

;; But it also looks as though clojure has been speeding up slightly, which has almost made up for that.

;; I'm also assuming that the JVM itself hasn't changed much since the original blog post.  The
;; netbook timings were all done today on the same JVM, but the desktop timings are from a while
;; back, so that might account for some differences.

;; So we're looking at a slowdown of about 300 times over what we could probably
;; achieve coding in assembler or in C with a good optimizing compiler (and of
;; course I'm assuming that floating point operations take one cycle each)

;; This is about the sort of speed that you'd expect from a dynamic language
;; without any optimization or type hinting.

;; In the original blog post I went through a number of faster versions:

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;; In solveit-2 I explicitly typed the loop variables, inlined the function call, and made an
;; internal target for recur. That speeded things up considerably:

(defn solveit-2 [t0 y0 h its]
  (loop [t0 (double t0), y0 (double y0), h (double h), its (int its)]
    (if (> its 0) 
      (let [t1 (+ t0 h)
            y1 (+ y0 (* h (- t0 y0)))]
        (recur t1 y1 h (dec its)))
      [t0 y0 h its])))

;; Let's time that and see how it goes:
(cyclesperit (solveit-2 0.0 1.0) 10000000) 

;; On my old desktop with clojure 1.2 : 490 cycles

;; On my netbook with clojure 1.2 : 600 cycles
;; On my netbook with clojure 1.3 : 44 cycles
;; On my netbook with clojure 1.4 : 44 cycles


;; Wow!

;; Again we see that my netbook is around 20% slower, but the performance change between clojure 1.2
;; and clojure 1.3 was incredible.

;; There's now about a factor of 6 between the clojure version and the imaginary assembler program
;; on the imaginary cpu in my head. And all I've done is to declare the types and inline the functions.

;; This is actually so impressive that I want to examine the three things in detail, to see how the
;; three things combine to cause the speedup. As I remember from clojure 1.2, you needed all three
;; changes to see any great difference. In 1.4 you can type hint separately from making an internal loop target.
;; All timings with clojure 1.4.0.

(defn f [t y] (- t y))

(defn solveit [t0 y0 h its]
  (if (> its 0) 
    (let [t1 (+ t0 h)
          y1 (+ y0 (* h (f t0 y0)))]
      (recur t1 y1 h (dec its)))
    [t0 y0 h its]))

(cyclesperit (solveit 0.0 1.0) 10000000) ; 2481

;; type hints: 

(defn ^double tf [^double t ^double y] (- t y))

(defn solveit-1-1 [^double t0 ^double y0 ^double h ^long its]
  (if (> its 0) 
    (let [t1 (+ t0 h)
          y1 (+ y0 (* h (tf t0 y0)))]
      (recur t1 y1 h (dec its)))
    [t0 y0 h its]))

(cyclesperit (solveit-1-0 0.0 1.0) 10000000) ; 289

;; Inline f
(defn solveit-1-2 [ ^double t0 ^double y0 ^double h ^long its]
  (if (> its 0) 
    (let [t1 (+ t0 h)
          y1 (+ y0 (* h (- t0 y0)))]
      (recur t1 y1 h (dec its)))
    [t0 y0 h its]))

(cyclesperit (solveit-1-2 0.0 1.0) 10000000) ; 44

;; Internal loop target
(defn solveit-1-3 [ ^double t0 ^double y0 ^double h ^long its]
  (loop [t0 t0, y0 y0, h h, its its]
    (if (> its 0) 
      (let [t1 (+ t0 h)
            y1 (+ y0 (* h (- t0 y0)))]
        (recur t1 y1 h (dec its)))
      [t0 y0 h its])))

(cyclesperit (solveit-1-3 0.0 1.0) 10000000) ; 44

;; Original version
(defn solveit-2 [t0 y0 h its]
  (loop [t0 (double t0), y0 (double y0), h (double h), its (int its)]
    (if (> its 0) 
      (let [t1 (+ t0 h)
            y1 (+ y0 (* h (- t0 y0)))]
        (recur t1 y1 h (dec its)))
      [t0 y0 h its])))

;; Let's time that and see how it goes:
(cyclesperit (solveit-2 0.0 1.0) 10000000) ; 44


;; This is pretty awesome. Type hints alone are giving us a factor of eight speedup, and inlining
;; the function then gives us another factor of 6. The internal loop target, which was originally a
;; bit of a hack to allow clojure to deduce types doesn't make any difference to the hinted version,
;; but in fact that still works, and gives the compiler enough information to get the same speed.

;; So far I really could not be more impressed.

;; In the original post, I then used a profiler to find out where the loop was running slow, and 
;; did some strange things to make it fast:

(defn solveit-3 [t0 y0 h its]
  (loop [t0 (double t0), y0 (double y0), h (double h), its (int its)]
    (if (> its (int 0)) 
      (let [t1 (+ t0 h)
            y1 (+ y0 (* h (- t0 y0)))]
        (recur t1 y1 h (dec its)))
      [t0 y0 h its])))

(cyclesperit (solveit-3 0.0 1.0) 10000000) 

;; On my old desktop with clojure 1.2 : 70 cycles

;; On my netbook with clojure 1.2 : 90 cycles
;; On my netbook with clojure 1.3 : 100 cycles
;; On my netbook with clojure 1.4 : 100 cycles

;; Originally, that (int 0) instead of 0 in the comparison made all the difference, causing a huge
;; and unexpected speedup.

;; The odd thing here is that the same program actually runs slightly slower in clojure 1.3/1.4 than
;; it does in clojure 1.2. Not much, and they're all considerably slower than clojure 1.4's best
;; shot so far.

;; What I think is happening here is that that int cast is taking up unnecessary time in the later clojures,
;; whereas in 1.2 it was the final piece of the puzzle as far as the compiler was concerned.

;; To test this:
(defn solveit-3-a [ t0 y0 h its]
  (loop [t0 (double t0), y0 (double y0), h (double h), its (long its)]
    (if (> its (long 0 )) 
      (let [t1 (+ t0 h)
            y1 (+ y0 (* h (- t0 y0)))]
        (recur t1 y1 h (dec its)))
      [t0 y0 h its])))

(cyclesperit (solveit-3-a 0.0 1.0) 10000000) ; 44

;; solveit-3-a runs at the same 44 cycles/iteration as we've been seeing before. I get the
;; impression that modern clojure prefers longs to ints.


;; The final craziness was to take the (int 0) out of the loop entirely:

(defn solveit-4 [t0 y0 h its]
  (let [zero (int 0)]
    (loop [t0 (double t0) y0 (double y0) h (double h) its (int its)]
      (if (> its zero) 
        (let [t1 (+ t0 h)
              y1 (+ y0 (* h (- t0 y0)))]
          (recur t1 y1 h (dec its)))
      [t0 y0 h its]))))

(cyclesperit (solveit-4 0.0  1.0) 100000000)

;; And this is where it gets weird:

;; On my old desktop with clojure 1.2 : 23 cycles

;; On my netbook with clojure 1.2 : 32 cycles
;; On my netbook with clojure 1.3 : 45 cycles
;; On my netbook with clojure 1.4 : 45 cycles

;; Clojure 1.3/1.4 is now running at its usual (fast) speed, despite the ints rather than longs in solveit-4

;; But Clojure 1.2 on this program is able to run about 33% faster than the modern versions.

;; So it actually looks as though whatever awesomeness has caused the phenomenal speedup between
;; clojure 1.2 and 1.3 has actually slightly slowed down the optimized version.

;; I'm imagining that this is something to do with using longs rather than ints. 

;; But I haven't done any kind of close profiling to see if I can make 1.4 run any faster than 44 cycles/loop.

;; Someone in a comment on my earlier post said that they were seeing 8 cycles/second on their
;; Macintosh, so it might just be that there's something weird about the Atom processor.

;; At any rate, it looks as though optimizing numeric code in Clojure is now dead easy. You just
;; tell it the types, inline function calls, and that's as good as it gets.

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;; For those in real need of cycles and willing to take risks to save them, there's:

(set! *unchecked-math* true)

;; Which I don't think has an equivalent in 1.2, which cuts a few cycles off:

(defn solveit-4 [t0 y0 h its]
  (let [zero (long 0)]
    (loop [t0 (double t0) y0 (double y0) h (double h) its (long its)]
      (if (> its zero) 
        (let [t1 (+ t0 h)
              y1 (+ y0 (* h (- t0 y0)))]
          (recur t1 y1 h (dec its)))
      [t0 y0 h its]))))

(cyclesperit (solveit-4 0.0  1.0) 100000000) ; 37

;; I do wonder what is going on here. I'm not sure what checking is here to be turned off.

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;; Incidentally, the following two java programs also produce the numbers that are seen above:

;; public class Euler{
;;     public static void main (String[] args){
;;         double cpuspeed = 1.662;
;;         int its = 10000000;
;;         double t=0;
;;         double y=0;
;;         double h=1.0/its;

;;         Long start = System.nanoTime();
;;         for(int i = 0; i < its; i++){
;;             y = y+h*(t-y);
;;             t = t+h;
;;         }
;;         Long finish = System.nanoTime();

;;         System.out.println("y=" + y + " t=" +t); 
;;         System.out.println("cycles/iteration: " +  ((int) ((cpuspeed * (finish - start)) / its)));
;;     }
;; }


;; y=0.3678794227282174 t=0.99999999975017
;; cycles/iteration: 32


;; public class Euler{
;;     public static void main (String[] args){
;;         double cpuspeed = 1.662;
;;         long its = 10000000;
;;         double t=0;
;;         double y=0;
;;         double h=1.0/its;

;;         long start = System.nanoTime();
;;         for(long i = 0; i < its; i++){
;;             y = y+h*(t-y);
;;             t = t+h;
;;         }
;;         long finish = System.nanoTime();

;;         System.out.println("y=" + y + " t=" +t); 
;;         System.out.println("cycles/iteration: " +  ((int) ((cpuspeed * (finish - start)) / its)));
;;     }
;; }

;; y=0.3678794227282174 t=0.99999999975017
;; cycles/iteration: 37

FOR vs MAP

;; List Comprehension vs Map and Filter

;; The consequences of map taking many possible arguments were not clear to me until I realized that
;; it's 'side by side' iteration.

;; I asked this question:
;; http://stackoverflow.com/questions/14737391/how-do-i-combine-two-vectors-of-vectors-element-wise-in-clojure
;; on Stack Overflow, and the answer made me feel like a complete idiot.

;; ---------------------------------------------------------------------------------------

;; I can write:

(for [i [1 2 3]] (* 2 i)) ; (2 4 6)

;; Or I can write 

(map (fn[i] (* 2 i)) [1 2 3]) ; (2 4 6)

;; And which is easier to read is a bit of a matter of taste

;; I can write 

(for [i [1 2 3]] (for  [j [4 5 6]] (* i j))) ; ((4 5 6) (8 10 12) (12 15 18))

;; Or I can write 

(map (fn[j] (map (fn[i] (* i j)) [4 5 6])) [1 2 3]) ; ((4 5 6) (8 10 12) (12 15 18))

;; And I think most people will find the first one easier to read. Certainly I found it easier to write.

;; I can write 
(for [i [1 2 3], j [4 5 6]] (* i j)) ; (4 5 6 8 10 12 12 15 18)

;; Or I can write 
(mapcat (fn[j] (map (fn[i] (* i j)) [4 5 6])) [1 2 3]) ; (4 5 6 8 10 12 12 15 18)

;; And so I had a sort of feeling that list comprehension and map/filter were equivalent but that
;; list comprehensions were somewhat easier to read and write.

;; This certainly seems to be the case in python, where there have been occasional attempts to take
;; map, filter and reduce out of the language, in favour of list comprehensions


;; But occasionally, I'd like to iterate over two structures at once.

;; So [1 2 3], [4 5 6], * -> [4 10 18], say

;; And I am usually writing some horror like:

(for [[i,j] (partition 2 (interleave [1 2 3] [4 5 6]))] (* i j)) ; (4 10 18)

;; But in fact:

(map * [1 2 3] [4 5 6]) ; (4 10 18)

;; In which case I am thinking that the version with map is way clearer, and I am also thinking that
;; perhaps list comprehensions cannot do everything that map and filter can.

;; Or maybe I just don't know the syntax for side by side iteration. Or even the proper name. The
;; earlier combination is called 'cartesian product', and I am tempted to call the latter one
;; 'elementwise combination'.

;; Can anyone enlighten me on this subject?

;; Here's a nice application of multi-argument map, closely related to a real problem I'm trying to solve:

(def a '[[c c c]
         [y y y]
         [m m m]])

(def b '[[r g b]
         [r g b]
         [r g b]])

(def c '[[x y z]
         [y z x]
         [z x y]])

(mapv (partial mapv vector) a b c) 
;  [[[c r x] [c g y] [c b z]] 
;   [[y r y] [y g z] [y b x]] 
;   [[m r z] [m g x] [m b y]]]

;; Thank you to Alex Stoddard for that, and indeed a certain amount of associated enlightenment!

Runtime Require: Downloading a New Dependency and Adding it to the Classpath at Runtime with Pomegranate

Everything that I don't like about Clojure is really something I don't like about the JVM. Here is one tiny blow for freedom:

;; How to pull in a library at runtime

;; dependencies [[org.clojure/clojure "1.4.0"]
;;               [com.cemerick/pomegranate "0.0.13"]]


;; When working at the REPL, it's often a pain to stop clojure, add a library to your dependencies, and restart.

;; Behold the solution, due to Chas Emerick (on whom blessings)

(require 'cemerick.pomegranate)

;; We don't, I think, even need to download the library. Pomegranate will go off and get it if it needs to:

(cemerick.pomegranate/add-dependencies 
 :coordinates '[[ring "1.1.7"]]
 :repositories {"clojars" "http://clojars.org/repo" } )

;; It takes a while, but once it's done :

(require 'ring.adapter.jetty)

(defonce server (ring.adapter.jetty/run-jetty (fn[x]{:body"<h1>yo</h1>"}) {:port 8080 :join? false}))

;; And it works!

;; I am told that there are reasons not to use this in production
;; code, only at the REPL, which is a terrible shame, since otherwise
;; programs could specify their own dependencies and versions rather
;; than having to exile them to a separate file.

;; Still, this is a huge time-saver, and takes a lot of the pain out of trying a new library.