Reaction-Diffusion Morphogenesis

canvasmathsimulation

In 1952, Alan Turing wrote a paper proposing a reaction–diffusion model as the basis of the development of patterns such as the spots and stripes seen in animal skin. Mathematically, reaction–diffusion systems take the form of semi-linear parabolic partial differential equations. By iteration 300, it should be clear how the two elements have separated from an initial randomized system.

19 Nov 2014 by rm-hull

rm-hull/07d0706f21551e4249ce

REACTION-DIFFUSION-MORPHOGENESIS.md

Different presets:

reaction-diffusion-morphogenesis.cljs

;; Turing's reaction-diffusion model.
;; For more information, see:
;;   [1] Rafael Collantes. Algorithm Alley. Dr. Dobb's Journal, December 1996.
;;   [2] Alan M. Turing. The chemical basis of morphogenesis. Philosophical
;;       Transactions of the Royal Society of London. B 327, 37–72 (1952)
;;   [3] http://www.cgjennings.ca/toybox/turingmorph/

(ns enchilada.reaction-diffusion-morphogenesis
  (:require
    [enchilada :refer [canvas ctx value-of canvas-size proxy-request]]
    [cljs.core.async :as async :refer [<!]]
    [jayq.core :refer [show]]
    [big-bang.core :refer [big-bang]]
    [inkspot.color :as c]
    [inkspot.color-chart :as cc]
    [monet.canvas :refer [fill-style text]]))

(defn gaussian-seq []
  ; See Knuth, ACP, Section 3.4.1 Algorithm C.
  (let [v1 (dec (* 2.0 (rand)))
        v2 (dec (* 2.0 (rand)))
        s  (+ (* v1 v1) (* v2 v2))]
    (if (or (>= s 1.0) (= s 0.0))
      (recur)
      (let [mult (Math/sqrt (* -2.0 (/ (Math/log s) s)))]
        (lazy-cat [(* v1 mult) (* v2 mult)] (gaussian-seq))))))

(defn randomize [width height]
  (->>
    (gaussian-seq)
    (take (* width height))
    (mapv #(+ (* (rand) 12.0) (* % 2.0)))))

(defn zeros [width height]
  (vec (repeat (* width height) 0)))

(def presets
  {:cheetah     {:ca 3.5 :cb 16.0}
   :colony      {:ca 1.6 :cb  6.0}
   :fine        {:ca 0.1 :cb  1.0}
   :fingerprint {:ca 1.0 :cb 16.0}
   :maze        {:ca 2.6 :cb 24.0}
   :pocked      {:ca 1.0 :cb  3.0}})

(def gradients
  [:dance-to-forget
   :behongo
   :dracula
   :frozen
   :misty-meadow
   :opa
   :shore
   :starfall
   :vasily
   :virgin
   :winter ])

(defn opaque [rgba]
  (assoc rgba 3 0xff))

(defn initial-state [width height]
  (color-mapper
    (solve
      (let [{:keys [ca cb]} (presets (keyword (value-of :preset :fingerprint)))]
        {:t 0
         :ca ca
         :cb cb
         :ao (randomize width height)
         :an (zeros width height)
         :bo (randomize width height)
         :bn (zeros width height)
         :width width
         :height height
         :gradient (mapv
                     (comp opaque c/string->color)
                     (cc/ui-gradient (rand-nth gradients) 256))
       }))))

(defn swap-buffers [{:keys [ao an bo bn] :as world-state}]
  (assoc world-state
    :ao an
    :an ao
    :bo bn
    :bn bo))

(defn offsets [i width height]
  (let [m (mod i width)
        h (* width (dec height))]
    [(if (< i width) (+ i h) (- i width))       ; NORTH
     (inc (if (= m (dec width)) (- i width) i)) ; EAST
     (if (>= i h) (- i h) (+ i width))          ; SOUTH
     (dec (if (zero? m) (+ i width) i))]))      ; WEST

(defn solve [{:keys [ca cb ao an bo bn width height] :as world-state}]
  (let [an (transient an)
        bn (transient bn)]
    (dotimes [i (* width height)]
      (let [[n e s w] (offsets i width height)
            di-a (* ca (+ (ao n) (ao e) (ao s) (ao w) (* -4.0 (ao i))))
            re-a (- (* (ao i) (bo i)) (ao i) 12.0)

            di-b (* cb (+ (bo n) (bo e) (bo s) (bo w) (* -4.0 (bo i))))
            re-b (- 16.0 (* (ao i) (bo i)))]
        (assoc! an i (max 0.0 (+ (ao i) (* 0.01 (+ re-a di-a)))))
        (assoc! bn i (max 0.0 (+ (bo i) (* 0.01 (+ re-b di-b)))))))
    (assoc world-state
      :an (persistent! an)
      :bn (persistent! bn))))

(defn minmax [[x & xs]]
  (reduce
    (fn [[tiny big] n] [(min tiny n) (max big n)])
    [x x]
    xs))

(defn color-mapper [{:keys [an gradient] :as world-state}]
  (assoc world-state
    :color-mapper (apply cc/color-mapper gradient (minmax an))))

(defn set-rgba! [pixels offset rgba]
  (when-not (empty? rgba)
    (aset pixels offset (first rgba))
    (recur pixels (inc offset) (next rgba))))

(defn render [{:keys [t width height an color-mapper]}]
  (let [img  (.getImageData ctx 0 0 width height)
        data (.-data img)]

    (dotimes [i (* width height)]
      (set-rgba!
        data
        (* i 4)
        (color-mapper (an i))))

    (dotimes [y 6]
      (dotimes [x 8]
        (.putImageData ctx img (* x 100) (* y 100))))

    (->
      ctx
      (fill-style :#333)
      (text {:text (str "Iteration: " t) :x 10 :y 10}))))

(defn update-state [event world-state]
  (-> world-state (update-in [:t] inc) swap-buffers solve color-mapper))

(show canvas)

(big-bang
  :initial-state (initial-state 100 100)
  :on-tick update-state
  :to-draw render)