New Gas Particle Sandbox

globals [
  tick-advance-amount               ; how much we advance the tick counter this time through
  max-tick-advance-amount           ; the largest tick-advance-amount is allowed to be
  init-avg-speed init-avg-energy    ; initial averages
  avg-speed avg-energy              ; current average
  avg-energy-green
  avg-energy-orange
  avg-energy-purple
  particle-size
  toggle-red-state
  toggle-green-state
  min-particle-energy
  max-particle-energy
  particles-to-add
  show-wall-hits?
  max-particles
  #particles
  #-green-particles
  #-purple-particles
  #-orange-particles
]

breed [ particles particle ]
breed [ walls wall ]
breed [ flashes flash ]
breed [ erasers eraser ]
breed [ arrowheads arrowhead ]

erasers-own [ pressure? ]
flashes-own [ birthday ]

particles-own [
  speed mass energy          ; particles info
  last-collision
  color-type
]

walls-own [
  energy
  valve-1?
  valve-2?
  pressure?
  surface-energy
]

to setup
  clear-all
  reset-ticks
  set particle-size 1.0
  set max-tick-advance-amount 0.02
  set show-wall-hits? false
  set particles-to-add 2

  set-default-shape flashes "square"
  set-default-shape walls "wall"
  set-default-shape erasers "eraser"
  set-default-shape arrowheads "default"
  set max-particles 200
  set #-green-particles 0
  set #-purple-particles 0
  set #-orange-particles 0

  set min-particle-energy 0
  set max-particle-energy 500  ;(.5 ) * ( max-dist-in-tick-advance-amount  / max-tick-advance-amount ) ^ 2

  create-erasers 1 [
    set hidden? true
    set pressure? true
    set size 3
    set color white
  ]

  make-box
  make-particles

  ask particles [ apply-speed-visualization ]

  set init-avg-speed avg-speed
  set init-avg-energy avg-energy

  update-variables
  do-plotting
end 

to go
  mouse-action
  if let-particles-move?   [
    ask particles [ bounce ]
    ask particles [ move ]
    ask particles [ check-for-collision ]
    ask particles with [ any? walls-here ] [ rewind-to-bounce ]
    ask particles with [ any? walls-here ] [ remove-from-walls ]
  ]

  tick-advance tick-advance-amount
  calculate-tick-advance-amount

  ask flashes [ apply-flash-visualization ]
  ask particles [ apply-speed-visualization ]

  update-variables
  do-plotting
  display
end 

to update-variables
  let green-particles particles with [ color-type = 55 ]
  let purple-particles  particles with [ color-type = 115 ]
  let orange-particles particles with [ color-type = 25 ]
  set #-green-particles count green-particles
  set #-orange-particles count orange-particles
  set #-purple-particles count purple-particles
  set #particles (#-green-particles + #-purple-particles + #-orange-particles)

  if any? particles [
    set avg-speed  mean [ speed ] of particles
    set avg-energy mean [ energy ] of particles
  ]

  ifelse #-green-particles > 0 [
    set avg-energy-green mean [ energy ] of green-particles
  ]
  [ set avg-energy-green "N/A"]

  ifelse  #-orange-particles > 0 [
    set avg-energy-orange mean [ energy ] of orange-particles
  ]
  [ set avg-energy-orange "N/A"]

  ifelse  #-purple-particles > 0 [
    set avg-energy-purple mean [ energy ] of purple-particles
  ]
  [ set avg-energy-purple "N/A"]
end 

to follow-a-particle
  if #particles > 0 [
    follow one-of particles
  ]
end 

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;WALL INTERACTION;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

to toggle-red-wall
  ifelse toggle-red-state = "closed" [
    ask walls with [ valve-1? ] [
      set hidden? true
    ]
    set toggle-red-state "open"
  ][
    ask walls with [ valve-1? ] [
      set hidden? false
    ]
    set toggle-red-state "closed"
  ]
end 

to toggle-green-wall
  ifelse toggle-green-state = "closed" [
    ask walls with [ valve-2? ] [
      set hidden? true
    ]
    set toggle-green-state "open"
  ][
    ask walls with [ valve-2?] [
      set hidden? false
    ]
    set toggle-green-state "closed"
  ]
end 

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;GAS MOLECULES MOVEMENT;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

to bounce  ; particles procedure
  ; get the coordinates of the patch we'll be on if we go forward 1
  let bounce-patch nobody
  let bounce-patches nobody
  let hit-angle 0
  let this-patch patch-here
  let new-px 0
  let new-py 0
  let visible-wall nobody

  set bounce-patch  min-one-of walls in-cone ((sqrt (2)) / 2) 180 with [ myself != this-patch ] [ distance myself ]

  if bounce-patch != nobody [
    set new-px [ pxcor ] of bounce-patch
    set new-py [ pycor ] of bounce-patch
    set visible-wall walls-on bounce-patch

    if any? visible-wall with [ not hidden? ]  [
      set hit-angle towards bounce-patch
      ifelse (hit-angle <= 135 and hit-angle >= 45) or (hit-angle <= 315 and hit-angle >= 225) [
        set heading (- heading)
      ][
        set heading (180 - heading)
      ]
      if show-wall-hits? [
        ask patch new-px new-py [
          sprout 1 [
            set breed flashes
            set color gray - 2
            set birthday ticks
          ]
        ]
      ]
    ]
  ]
end 

to rewind-to-bounce  ; particles procedure
  ; attempts to deal with particle penetration by rewinding the particle path back to a point
  ; where it is about to hit a wall
  ; the particle path is reversed 49% of the previous tick-advance-amount it made,
  ; then particle collision with the wall is detected again.
  ; and the particle bounces off the wall using the remaining 51% of the tick-advance-amount.
  ; this use of slightly more of the tick-advance-amount for forward motion off the wall, helps
  ; insure the particle doesn't get stuck inside the wall on the bounce.

  let bounce-patch nobody
  let bounce-patches nobody
  let hit-angle 0
  let this-patch nobody
  let new-px 0
  let new-py 0
  let visible-wall nobody

  bk (speed) * tick-advance-amount * .49
  set this-patch  patch-here

  set bounce-patch  min-one-of walls in-cone ((sqrt (2)) / 2) 180 with [ self != this-patch ] [ distance myself ]

  if bounce-patch != nobody [

    set new-px [pxcor] of bounce-patch
    set new-py [pycor] of bounce-patch
    set visible-wall walls-on bounce-patch

    if any? visible-wall with [not hidden?] [
      set hit-angle towards bounce-patch

      ifelse (hit-angle <= 135 and hit-angle >= 45) or (hit-angle <= 315 and hit-angle >= 225) [
        set heading (- heading)
      ][
        set heading (180 - heading)
      ]

      if show-wall-hits? [
        ask patch new-px new-py [
          sprout 1 [
            set breed flashes
            set color gray - 2
            set birthday ticks
          ]
        ]
      ]
    ]
  ]
  fd (speed) * tick-advance-amount * 2
end 

to move  ; particles procedure
  if patch-ahead (speed * tick-advance-amount) != patch-here [ set last-collision nobody ]
  fd (speed * tick-advance-amount * 2)
end 

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;GAS MOLECULES COLLISIONS;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;from GasLab

to calculate-tick-advance-amount
  ; tick-advance-amount is calculated in such way that even the fastest
  ; particles will jump at most 1 patch delta in a ticks tick. As
  ; particles jump (speed * tick-advance-amount) at every ticks tick, making
  ; tick delta the inverse of the speed of the fastest particles
  ; (1/max speed) assures that. Having each particles advance at most
  ; one patch-delta is necessary for it not to "jump over" a wall
  ; or another particles.
  ifelse any? particles with [ speed > 0 ] [
    set tick-advance-amount min list (1 / (ceiling max [speed] of particles )) max-tick-advance-amount
  ][
    set tick-advance-amount max-tick-advance-amount
  ]
end 

to check-for-collision  ; particles procedure
  ; Here we impose a rule that collisions only take place when there
  ; are exactly two particles per patch.  We do this because when the
  ; student introduces new particles from the side, we want them to
  ; form a uniform wavefront.
  ;
  ; Why do we want a uniform wavefront?  Because it is actually more
  ; realistic.  (And also because the curriculum uses the uniform
  ; wavefront to help teach the relationship between particles collisions,
  ; wall hits, and pressure.)
  ;
  ; Why is it realistic to assume a uniform wavefront?  Because in reality,
  ; whether a collision takes place would depend on the actual headings
  ; of the particles, not merely on their proximity.  Since the particles
  ; in the wavefront have identical speeds and near-identical headings,
  ; in reality they would not collide.  So even though the two-particles
  ; rule is not itself realistic, it produces a realistic result.  Also,
  ; unless the number of particles is extremely large, it is very rare
  ; for three or  particles to land on the same patch (for example,
  ; with 400 particles it happens less than 1% of the time).  So imposing
  ; this additional rule should have only a negligible effect on the
  ; aggregate behavior of the system.
  ;
  ; Why does this rule produce a uniform wavefront?  The particles all
  ; start out on the same patch, which means that without the only-two
  ; rule, they would all start colliding with each other immediately,
  ; resulting in much random variation of speeds and headings.  With
  ; the only-two rule, they are prevented from colliding with each other
  ; until they have spread out a lot.  (And in fact, if you observe
  ; the wavefront closely, you will see that it is not completely smooth,
  ; because  collisions eventually do start occurring when it thins out while fanning.)

  if count other particles-here  in-radius 1 = 1 [
    ; the following conditions are imposed on collision candidates:
    ;   1. they must have a lower who number than my own, because collision
    ;      code is asymmetrical: it must always happen from the point of view
    ;      of just one particles.
    ;   2. they must not be the same particles that we last collided with on
    ;      this patch, so that we have a chance to leave the patch after we've
    ;      collided with someone.
    let candidate one-of other particles-here with [ who < [ who ] of myself and myself != last-collision ]
    ;; we also only collide if one of us has non-zero speed. It's useless
    ;; (and incorrect, actually) for two particles with zero speed to collide.
    if (candidate != nobody) and (speed > 0 or [ speed ] of candidate > 0) [
      collide-with candidate
      set last-collision candidate
      ask candidate [ set last-collision myself ]
    ]
  ]
end 

; implements a collision with another particles.
;
; THIS IS THE HEART OF THE particles SIMULATION, AND YOU ARE STRONGLY ADVISED
; NOT TO CHANGE IT UNLESS YOU REALLY UNDERSTAND WHAT YOU'RE DOING!
;
; The two particles colliding are self and other-particles, and while the
; collision is performed from the point of view of self, both particles are
; modified to reflect its effects. This is somewhat complicated, so I'll
; give a general outline here:
;   1. Do initial setup, and determine the heading between particles centers
;      (call it theta).
;   2. Convert the representation of the velocity of each particles from
;      speed/heading to a theta-based vector whose first component is the
;      particle's speed along theta, and whose second component is the speed
;      perpendicular to theta.
;   3. Modify the velocity vectors to reflect the effects of the collision.
;      This involves:
;        a. computing the velocity of the center of mass of the whole system
;           along direction theta
;        b. updating the along-theta components of the two velocity vectors.
;   4. Convert from the theta-based vector representation of velocity back to
;      the usual speed/heading representation for each particles.
;   5. Perform final cleanup and update derived quantities.

to collide-with [ other-particles ] ;; particles procedure
  ; PHASE 1: initial setup

  ; for convenience, grab  quantities from other-particles
  let mass2 [ mass ] of other-particles
  let speed2 [ speed ] of other-particles
  let heading2 [ heading ] of other-particles

  ; since particles are modeled as zero-size points, theta isn't meaningfully
  ; defined. we can assign it randomly without affecting the model's outcome.
  let theta (random-float 360)

  ; PHASE 2: convert velocities to theta-based vector representation

  ; now convert my velocity from speed/heading representation to components
  ; along theta and perpendicular to theta
  let v1t (speed * cos (theta - heading))
  let v1l (speed * sin (theta - heading))

  ;; do the same for other-particles
  let v2t (speed2 * cos (theta - heading2))
  let v2l (speed2 * sin (theta - heading2))

  ; PHASE 3: manipulate vectors to implement collision

  ; compute the velocity of the system's center of mass along theta
  let vcm (((mass * v1t) + (mass2 * v2t)) / (mass + mass2) )

  ; now compute the new velocity for each particles along direction theta.
  ; velocity perpendicular to theta is unaffected by a collision along theta,
  ; so the next two lines actually implement the collision itself, in the
  ; sense that the effects of the collision are exactly the following changes
  ; in particles velocity.
  set v1t (2 * vcm - v1t)
  set v2t (2 * vcm - v2t)

  ; PHASE 4: convert back to normal speed/heading

  ; now convert my velocity vector into my new speed and heading
  set speed sqrt ((v1t ^ 2) + (v1l ^ 2))
  set energy (0.5 * mass * speed ^ 2)
  ; if the magnitude of the velocity vector is 0, atan is undefined. but
  ; speed will be 0, so heading is irrelevant anyway. therefore, in that
  ; case we'll just leave it unmodified.
  if v1l != 0 or v1t != 0 [ set heading (theta - (atan v1l v1t)) ]

  ;; and do the same for other-particle
  ask other-particles [
    set speed sqrt ((v2t ^ 2) + (v2l ^ 2))
    set energy (0.5 * mass * (speed ^ 2))
    if v2l != 0 or v2t != 0 [ set heading (theta - (atan v2l v2t)) ]
  ]
end 



;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;  mouse interaction procedures
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

to mouse-action
 let snap-xcor 0
 let snap-ycor 0
 let orig-xcor 0
 let orig-ycor 0
 let eraser-window-walls nobody
 let eraser-window-particles nobody

  ifelse mouse-down? [
    set orig-xcor mouse-xcor
    set orig-ycor mouse-ycor
    set snap-xcor round orig-xcor
    set snap-ycor round orig-ycor

    ask patches with [ pxcor = snap-xcor and pycor = snap-ycor ] [
      set eraser-window-walls (turtle-set walls-here walls-on neighbors)
      set eraser-window-walls eraser-window-walls with [ not pressure? ]
      set eraser-window-particles (turtle-set particles-here particles-on neighbors)

      if mouse-interaction = "draw basic wall" [
        ask walls-here [ die ]
        sprout 1 [
          set breed walls set color gray
          initialize-this-wall
        ]
      ]

      if mouse-interaction = "draw red removable wall"  [
        set toggle-red-state "open"
        toggle-red-wall
        ask walls-here [ die ]
        sprout 1 [
          set breed walls set color red
          initialize-this-wall set valve-1? true
        ]
      ]

      if mouse-interaction = "draw green removable wall" [
        set toggle-green-state "open"
        toggle-green-wall
        ask walls-here [ die ]
        sprout 1 [
          set breed walls set color 55
          initialize-this-wall set valve-2? true
        ]
      ]

      if mouse-interaction = "big eraser" [
        ask erasers [
          set hidden? false
          set shape "eraser"
          setxy orig-xcor orig-ycor
        ]
        ask eraser-window-walls [ die ]
        ask eraser-window-particles [ die ]
      ]

      if mouse-interaction = "add purple particles" or mouse-interaction = "add green particles"
        or mouse-interaction = "add orange particles" and #particles < max-particles [

        ask erasers [
          set hidden? false
          set shape "spray paint"
          setxy orig-xcor orig-ycor
        ]
        sprout particles-to-add [
          set breed particles
          setup-particles
          jump random-float 2
          if mouse-interaction = "add purple particles" [ set color-type 115 color-particle-and-link ]
          if mouse-interaction = "add orange particles" [ set color-type 25 color-particle-and-link ]
          if mouse-interaction = "add green particles" [ set color-type 55  color-particle-and-link ]
          apply-speed-visualization
        ]
      ]

      if mouse-interaction = "paint particles purple"
        or mouse-interaction = "paint particles orange"
        or mouse-interaction = "paint particles green" [

        ask erasers [
          set hidden? false
          set shape "spray paint"
          setxy orig-xcor orig-ycor
        ]
        ask eraser-window-particles [
          if mouse-interaction = "paint particles purple"  [ set color-type 115 color-particle-and-link ]
          if mouse-interaction = "paint particles orange" [ set color-type 25 color-particle-and-link ]
          if mouse-interaction = "paint particles green" [ set  color-type 55 color-particle-and-link ]
          apply-speed-visualization
        ]
      ]

      if mouse-interaction = "speed up particles"  [

        ask erasers [
          set hidden? false
          set shape "spray paint"
          setxy orig-xcor orig-ycor
        ]
        ask eraser-window-particles [
          set energy (energy * 1.1)
          set energy limited-particle-energy
          set speed speed-from-energy
          apply-speed-visualization
        ]
      ]

      if mouse-interaction = "slow down particles" [
        ask erasers [
          set hidden? false
          set shape "spray paint"
          setxy orig-xcor orig-ycor
        ]
        ask eraser-window-particles [
          set energy (energy / 1.1)
          set energy limited-particle-energy
          set speed speed-from-energy
          apply-speed-visualization
        ]
      ]
    ]
    ask particles with [ any? walls-here ] [ remove-from-walls ] ; deal with any walls drawn on top of particles
  ][
    ask erasers [ set hidden? true ]
  ]
end 

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; particle speed and flash visualization procedures
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

to apply-flash-visualization
  set size (size * .95)
  if (ticks - birthday > 0.4)  [ die ]
end 

to apply-speed-visualization
  if visualize-speed? = "arrows" [ scale-arrowheads ]
  if visualize-speed? = "different shades" [ recolorshade ]
  if visualize-speed? = "none" [ recolornone ]
end 

to color-particle-and-link
  let this-link my-out-links
  let this-color-type color-type
  set color this-color-type
  ask this-link [ set color this-color-type ]
end 

to scale-arrowheads
  let this-xcor xcor
  let this-ycor ycor
  let this-speed speed
  let this-heading heading
  let this-arrowhead out-link-neighbors
  let this-link my-out-links
  ask this-link [ set hidden? false ]
  ask this-arrowhead [
    set xcor this-xcor
    set ycor this-ycor
    set heading this-heading
    fd .5 + this-speed * .5
  ]
end 

to recolorshade
  let this-link my-out-links
  ask this-link [ set hidden? true ]
  ifelse speed < 21 [
    set color color-type - 3 + speed  / 3
  ][
    set color color-type + 4
  ]
end 

to recolornone
  let this-link my-out-links
  ask this-link [ set hidden? true ]
  set color color-type
end 

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;  initialization procedures
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

to make-box
  ask patches with [ (pycor = min-pycor or pycor = max-pycor or pxcor = min-pxcor or pxcor = max-pxcor) ] [
    sprout 1 [
      set breed walls set color yellow initialize-this-wall set pressure? true
    ]
  ]
end 

to initialize-this-wall
  set valve-1? false
  set valve-2? false
  set pressure? false
end 

to make-particles
  create-particles initial-#-particles [
    setup-particles
    random-position
  ]
end 

to setup-particles  ; particles procedure
  set shape "circle"
  set size particle-size
  set energy initial-temperature
  set color-type 115
  set color color-type
  set mass (10)  ; atomic masses of oxygen atoms
  hatch 1 [
    set breed arrowheads
    set hidden? true
    create-link-from myself [ tie ]
  ]
  set speed speed-from-energy
  set last-collision nobody
end 


; Place particles at random, but they must not be placed on top of wall atoms.
; This procedure takes into account the fact that wall molecules could have two possible arrangements,
; i.e. high-surface area ot low-surface area.

to random-position ;; particles procedure
  let open-patches nobody
  let open-patch nobody
  set open-patches patches with [not any? turtles-here and pxcor != max-pxcor and pxcor != min-pxcor and pycor != min-pycor and pycor != max-pycor]
  set open-patch one-of open-patches

  ; Reuven added the following "if" so that we can get through setup without a runtime error.
  if open-patch = nobody [
    user-message "No open patches found.  Exiting."
    stop
  ]

  setxy ([ pxcor ] of open-patch) ([ pycor ] of open-patch)
end 


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; wall penetration error handling procedure
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; if particles actually end up within the wall

to remove-from-walls
  let this-wall walls-here with [ not hidden? ]

  if count this-wall != 0 [
    let available-patches patches with [ not any? walls-here ]
    let closest-patch nobody
    if (any? available-patches) [
      set closest-patch min-one-of available-patches [ distance myself ]
      set heading towards closest-patch
      setxy ([ pxcor ] of closest-patch)  ([ pycor ] of closest-patch)
    ]
  ]
end 


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;GRAPHS;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

to do-plotting

  set-current-plot "Temperature of gases"

  if any? particles with [ color-type = 55 ] [
    set-current-plot-pen "green"
    plotxy ticks avg-energy-green
  ]
  if any? particles with [ color-type = 25 ] [
    set-current-plot-pen "orange"
    plotxy ticks avg-energy-orange
  ]
  if any? particles with [color-type = 115] [
    set-current-plot-pen "purple"
    plotxy ticks avg-energy-purple
  ]
end 


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;REPORTERS;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

to-report speed-from-energy
  report sqrt (2 * energy / mass)
end 

to-report energy-from-speed
  report (mass * speed * speed / 2)
end 

to-report limited-particle-energy
  let limited-energy energy
  if limited-energy > max-particle-energy [ set limited-energy max-particle-energy ]
  if limited-energy < min-particle-energy [ set limited-energy min-particle-energy ]
  report limited-energy
end 


; Copyright 2006 Uri Wilensky.
; See Info tab for full copyright and license.