examples/interactions.rkt
#lang racket/gui
;;; interactions.rkt - Interactions Demo

(require (planet williams/simulation/simulation))
(require (planet williams/science/math))

;;; Simulation Parameters

;;; F : real?
;;; p-life : (>/c 0.0)
;;; n-min : exact-positive-integer?
;;; n-max : exact-positive-integer?
;;; x-min : real?
;;; x-max : real?
;;; y-min : real?
;;; y-max : real?
;;; v-min : real?
;;; v-max : real?
(define F 100.0)                        ; force constant
(define p-life 5000.0)                  ; maximum particle lifetime
(define-values (n-min n-max) (values 20 40))
(define-values (x-min x-max) (values -1000.0 1000.0))
(define-values (y-min y-max) (values -1000.0 1000.0))
(define-values (v-min v-max) (values -1.0 1.0))

;;; Global Variables - set in run-simulation

;;; n : exact-nonnegative-integer?
;;; k : (or/c variable? false/c)
;;; p : (or/c (vectorof process?) false/c)
;;; x : (or/c (vectorof real?) false/c)
;;; y : (or/c (vectorof real?) false/c)
;;; dx/dt : (or/c (vectorof real?) false/c)
;;; dy/dt : (or/c (vectorof real?) false/c)
(define n 0)                            ; total number of particles
(define k #f)                           ; simulation progress variable
(define p #f)                           ; particle vector
(define x #f)                           ; particle position x vector
(define y #f)                           ; particle position y vector
(define dx/dt #f)                       ; particle position dx/dt vector
(define dy/dt #f)                       ; particle position dy/dt vector
(define pen #f)                         ; pen to draw the particle

;;; Particle Model

;;; process (particle i) -> void?
;;;   i : exact-nonnegative-integer?
(define-process (particle i)
  (let ((death-time (+ (current-simulation-time) (random-flat 0.0 p-life))))
    ;; let others know we exist
    (vector-set! p i self)
    ;; initialize velocities
    (vector-set! dx/dt i (make-continuous-variable (random-flat v-min v-max)))
    (vector-set! dy/dt i (make-continuous-variable (random-flat v-min v-max)))
    ;; initialize positions
    (vector-set! x i (make-continuous-variable (random-flat x-min x-max)))
    (vector-set! y i (make-continuous-variable (random-flat y-min y-max)))
    ;; continuous model
    (work/continuously
      until (>= (current-simulation-time) death-time)
      (let ((xi (variable-value (vector-ref x i)))
            (yi (variable-value (vector-ref y i)))
            (x-dot 0.0)
            (y-dot 0.0))
        ;; compute accelerations
        (for ((ii (in-range n)))
          (when (and (not (= ii i))
                     (vector-ref p ii))
            (let* ((xii (variable-value (vector-ref x ii)))
                   (yii (variable-value (vector-ref y ii)))
                   (r2 (max 100.0 
                            (+ (sqr (- xii xi)) (sqr (- yii yi))))))
              (set! x-dot (+ x-dot (/ (* (sign (- xii xi)) F) r2)))
              (set! y-dot (+ y-dot (/ (* (sign (- yii yi)) F) r2))))))
        ;; set velocity derivatives (i.e. accelerations)
        (set-variable-dt! (vector-ref dx/dt i) x-dot)
        (set-variable-dt! (vector-ref dy/dt i) y-dot)
        ;; set position derivatives (i.e. velocities)
        (set-variable-dt! (vector-ref x i) (variable-value (vector-ref dx/dt i)))
        (set-variable-dt! (vector-ref y i) (variable-value (vector-ref dy/dt i)))))
    ;; let others know we no longer exist
    (vector-set! p i #f)
    (set-variable-value! k (+ (variable-value k) 1))))

;;; Simulation Model

;;; (run-simulation) -> void?
(define (run-simulation)
  (break-enabled #t)
  (with-new-simulation-environment
    ;; initialize graphics
   (begin-busy-cursor)
   (send run-button enable #f)
   (send (send canvas get-dc) clear)
    ;; randomize
    (random-source-randomize! (current-random-source))
    ;; fixed step size
    (current-simulation-step-size 1.0)
    (current-simulation-control #f)
    ;; set number of particles randomly
    (set! n (+ n-min (random-uniform-int (- n-max n-min))))
    (send gauge set-range n)
    (send gauge set-value 0)
    ;; create and monitor simulation progress variable, k
    (set! k (make-variable 0))
    (monitor after (set-variable-value! k v)
             (send gauge set-value v))
    ;; create global variable vectors
    (set! p (make-vector n #f))
    (set! x (make-vector n))
    (set! y (make-vector n))
    (set! dx/dt (make-vector n))
    (set! dy/dt (make-vector n))
    (set! pen (build-vector
               n
               (lambda (i)
                 (new pen%
                      (color (make-object color%
                               (random 256) (random 256) (random 256)))))))
    ;; create the particles
    (for ((i (in-range n)))
      (schedule #:now (particle i)))
    ;; set simulation monitor
    (current-simulation-monitor
     (lambda ()
       (when (check-for-break)
         (stop-simulation))
       (let ((dc (send canvas get-dc)))
         (for ((i (in-range n)))
           (when (vector-ref p i)
             (let-values (((sx sy) (model-coords->screen-coords
                                    (variable-value (vector-ref x i))
                                    (variable-value (vector-ref y i)))))
               (send dc set-pen (vector-ref pen i))
               (send dc draw-point sx sy)))))
       (yield)))
  ;; start the simulation
  (start-simulation)
  ;; finalize graphics
  (send run-button enable #t)
  (end-busy-cursor)))

;;; Simulation Graphics

;;; WIDTH : exact-positive-integer?
;;; HEIGHT : exact-positive-integer?
(define-values (WIDTH HEIGHT) (values 800 600))

;;; (model-coords->screen-coords x y) -> real? real?
;;;   x : real?
;;;   y : real?
;;; convert model coordinates to screen coordinates
(define (model-coords->screen-coords x y)
  (values (* WIDTH (/ (- x x-min) (- x-max x-min)))
          (* HEIGHT (/ (- y y-min) (- y-max y-min)))))

;;; Graphical Elements

;;; frame : (is-a/c frame%)
(define frame
  (instantiate frame%
    ("Interacting Particles")
    (x 0)
    (y 0)))

;;; cancas : (is-a/c canvas%)
(define canvas
  (instantiate canvas%
    (frame)
    (min-width WIDTH)
    (min-height HEIGHT)))

;;; gauge : (is-a/c gauge%)
(define gauge
  (instantiate gauge%
    ("Progress" n-max frame)))

;;; run-button : (is-a/c button?)
(define run-button
  (instantiate button%
    ("Run" frame
           (lambda (b e)
             (run-simulation)))))

;;; Show the frame.
(send frame show #t)