core.clj

(ns analysis.core)

(use 'criterium.core)

;; algorithm: see http://en.wikipedia.org/wiki/Hilbert_curve

;;(set! *unchecked-math* true)

(deftype Point [^long x ^long y])

(defn hilbert-to-point [^long order ^long d]
  (binding [*unchecked-math* true]
    (let [n (bit-shift-left 1 (long order))]
      (loop [s 1, rx 0, ry 0, t d, x 0, y 0]
        (let [rx (bit-and 1 (bit-shift-right t 1))
              ry (bit-and 1 (bit-xor t rx))
              sm (dec s)
              x (long (if (== 1 ry) x (if (== 1 rx) (- sm y) y)))
              y (long (if (== 1 ry) y (if (== 1 rx) (- sm x) x)))
              x (+ x (* s rx))
              y (+ y (* s ry))
              t (bit-shift-right t 2)
              s (bit-shift-left s 1)]
          (if (< s n) (recur s rx ry t x y)
              (Point. x y))))))
  )

(defn red [^long o]
  (binding [*unchecked-math* true]
    (let [n2 (bit-shift-left 1 (* 2 o))]
      (loop [sum (Point. 0 0)
             d (long 0)]
        (let [^Point a (hilbert-to-point o d)]
          (if (< d n2)
            (recur (Point. (+ (.x a) (.x sum))
                           (+ (.y a) (.y sum)))
                   (inc d) )
            sum))))))
  
(prn (red 10))

(bench (red 10) :verbose)