GasLab Circular Particles

GasLab Circular Particles preview image

1 collaborator

Uri_dolphin3 Uri Wilensky (Author)


(This model has yet to be categorized with any tags)
Model group CCL | Visible to everyone | Changeable by group members (CCL)
Model was written in NetLogo 5.0.4 • Viewed 124 times • Downloaded 10 times • Run 2 times
Download the 'GasLab Circular Particles' modelDownload this modelEmbed this model

Do you have questions or comments about this model? Ask them here! (You'll first need to log in.)


This model is one in a series of GasLab models. They use the same basic rules for simulating the behavior of gases. Each model integrates different features in order to highlight different aspects of gas behavior.

This model is different from the other GasLab models in that the collision calculations take the circular shape and size of the particles into account, instead of modeling the particles as dimensionless points.


The model determines the resulting motion of particles that collide, with no loss in their total momentum or total kinetic energy (an elastic collision).

To calculate the outcome of collision, it is necessary to calculate the exact time at which the edge of one particle (represented as a circle), would touch the edge of another particle (or the walls of a container) if the particles were allowed to continue with their current headings and speeds.

By performing such a calculation, one can determine when the next collision anywhere in the system would occur in time. From this determination, the model then advances the motion of all the particles using their current headings and speeds that far in time until this next collision point is reached. Exchange of kinetic energy and momentum between the two particles, according to conservation of kinetic energy and conservation of momentum along the collision axis (a line drawn between the centers of the two particles), is then calculated, and the particles are given new headings and speeds based on this outcome.


INITIAL-NUMBER-PARTICLES determines the number of gas particles used with SETUP. If the world is too small or the particles are too large, the SETUP procedure of the particles will stop so as to prevent overlapping particles.

SMALLEST-PARTICLE-SIZE and LARGEST-PARTICLE-SIZE determines the range of particle sizes that will be created when SETUP is pressed. (Particles are also assigned a mass proportional to the area of the particle that is created.)

The SETUP button will set the initial conditions.
The GO button will run the simulation.


  • % FAST, % MEDIUM, % SLOW: the percentage of particles with different speeds: fast (red), medium (green), and slow (blue).
  • AVERAGE SPEED: average speed of the particles.
  • AVERAGE ENERGY: average kinetic energy of the particles.


  • SPEED COUNTS: plots the number of particles in each range of speed.
  • SPEED HISTOGRAM: speed distribution of all the particles. The gray line is the average value, and the black line is the initial average.
  • ENERGY HISTOGRAM: distribution of energies of all the particles, calculated as m*(v^2)/2. The gray line is the average value, and the black line is the initial average.

Initially, all the particles have the same speed but random directions. Therefore the first histogram plots of speed will show only one column. If all the particles have the same size (and therefore the same mass), then the first histogram plot of energy will also show one column. As the particles repeatedly collide, they exchange energy and head off in new directions, and the speeds are dispersed -- some particles get faster, some get slower. The histogram distribution changes accordingly.


Run the model with different masses, by setting the MAX-PARTICLE-SIZE larger than the MIN-PARTICLE-SIZE. Mass is scaled linearly to the area of the particle, so that a particle that is twice the radius of another particle, has four the area and therefore four times the mass.

With many different mass particles colliding over time, different sized particles start to move at different speed ranges (in general). The smallest mass particles will be usually moving faster (red) than the average particle speed and the largest mass particles will be usually slower (blue) than the average particle speed. This emergent result is what happens in a gas that is a mixture of particles of different masses. At any given temperature, the higher mass particles are moving slower (such as Nitrogen gas: N2) then the lower mass particles (such as water vapor: H2O).

The particle histograms quickly converge on the classic Maxwell-Boltzmann distribution. What's special about these curves? Why is the shape of the energy curve not the same as the speed curve?

Look at the other GasLab models to compare the particle histograms. In those models, the particles are modeled as "point" particles, with an area of "zero". The results of those models should match this model, even with that simplification.

With particles of different sizes, you may notice some fast moving particles have lower energy than medium speed particles. How can the difference in the mass of the particles account for this?


Setting all the particles to have a very slow speed (e.g. 0.001) and one particle to have a very fast speed helps show how kinetic energy is eventually transferred to all the particles through a series of collisions and would serve as a good model for energy exchange through conduction between hot and cold gases.

To see what the approximate mass of each particle is, type this in the command center:

ask particles [ set label precision mass 0 ]


Collisions between boxes and circles could also be explored. Variations in size between particles could investigated or variations in the mass of some of the particle could be made to explore other factors that affect the outcome of collisions.


Instead of advancing one tick at a time as in most models, the tick counter takes on fractional values, using the tick-advance primitive. (In the Interface tab, it is displayed as an integer, but if you make a monitor for ticks you'll see the exact value.)


Look at the other GasLab models to see collisions of "point" particles, that is, the particles are assumed to have an area or volume of zero.


This model was developed as part of the GasLab curriculum ( and has also been incorporated into the Connected Chemistry curriculum (


If you mention this model in a publication, we ask that you include these citations for the model itself and for the NetLogo software:


Copyright 2005 Uri Wilensky.


This work is licensed under the Creative Commons Attribution-NonCommercial-ShareAlike 3.0 License. To view a copy of this license, visit or send a letter to Creative Commons, 559 Nathan Abbott Way, Stanford, California 94305, USA.

Commercial licenses are also available. To inquire about commercial licenses, please contact Uri Wilensky at

This model and associated activities and materials were created as part of the project: MODELING ACROSS THE CURRICULUM. The project gratefully acknowledges the support of the National Science Foundation, the National Institute of Health, and the Department of Education (IERI program) -- grant number REC #0115699.

Comments and Questions

Click to Run Model

globals [
  tick-delta        ;; how much simulation time will pass in this step
  box-edge          ;; distance of box edge from origin
  collisions        ;; list used to keep track of future collisions
  particle1         ;; first particle currently colliding
  particle2         ;; second particle currently colliding
  init-avg-speed init-avg-energy             ;; initial averages
  avg-speed avg-energy                       ;; current averages
  fast medium slow                           ;; current counts
  percent-slow percent-medium percent-fast   ;; percentage of current counts

breed [particles particle]

particles-own [

;;;;;;;;;;;;;;;;;; setup procedures ;;;;;;;;;;;;;;;;;;;;;

to setup
  set-default-shape particles "circle"
  set box-edge max-pxcor - 1
  ask patches with [(abs pxcor = box-edge or abs pycor = box-edge) and
                    abs pxcor <= box-edge and abs pycor <= box-edge]
    [ set pcolor yellow ]
  set avg-speed 1
  set particle1 nobody
  set particle2 nobody
  set collisions []
  ask particles [ check-for-wall-collision ]
  ask particles [ check-for-particle-collision ]
  set init-avg-speed avg-speed
  set init-avg-energy avg-energy

to make-particles
  create-particles initial-number-particles [
    set speed 1

    set size smallest-particle-size
             + random-float (largest-particle-size - smallest-particle-size)
    ;; set the mass proportional to the area of the particle
    set mass (size * size)
    set energy kinetic-energy

  ;; When space is tight, placing the big particles first improves
  ;; our chances of eventually finding places for all of them.
  foreach sort-by [[size] of ?1 > [size] of ?2] particles [
    ask ? [
      while [overlapping?] [ position-randomly ]

to-report overlapping?  ;; particle procedure
  ;; here, we use IN-RADIUS just for improved speed; the real testing
  ;; is done by DISTANCE
  report any? other particles in-radius ((size + largest-particle-size) / 2)
                              with [distance myself < (size + [size] of myself) / 2]

to position-randomly  ;; particle procedure
  ;; place particle at random location inside the box
  setxy one-of [1 -1] * random-float (box-edge - 0.5 - size / 2)
        one-of [1 -1] * random-float (box-edge - 0.5 - size / 2)

;;;;;;;;;;;;;;;;;;;;;;;; go procedures  ;;;;;;;;;;;;;;;;;

to go
  ask particles [ jump speed * tick-delta ]
  tick-advance tick-delta
  if floor ticks > floor (ticks - tick-delta)

to update-variables
  set medium count particles with [color = green]
  set slow   count particles with [color = blue]
  set fast   count particles with [color = red]
  set percent-medium (medium / ( count particles )) * 100
  set percent-slow (slow / (count particles)) * 100
  set percent-fast (fast / (count particles)) * 100
  set avg-speed  mean [speed] of particles
  set avg-energy mean [energy] of particles

to recalculate-particles-that-just-collided
  ;; Since only collisions involving the particles that collided most recently could be affected,
  ;; we filter those out of collisions.  Then we recalculate all possible collisions for
  ;; the particles that collided last.  The ifelse statement is necessary because
  ;; particle2 can be either a particle or a string representing a wall.  If it is a
  ;; wall, we don't want to invalidate all collisions involving that wall (because the wall's
  ;; position wasn't affected, those collisions are still valid.
  ifelse is-turtle? particle2
      set collisions filter [item 1 ? != particle1 and
                             item 2 ? != particle1 and
                             item 1 ? != particle2 and
                             item 2 ? != particle2]
      ask particle2 [ check-for-wall-collision ]
      ask particle2 [ check-for-particle-collision ]
      set collisions filter [item 1 ? != particle1 and
                             item 2 ? != particle1]
  if particle1 != nobody [ ask particle1 [ check-for-wall-collision ] ]
  if particle1 != nobody [ ask particle1 [ check-for-particle-collision ] ]
  ;; Slight errors in floating point math can cause a collision that just
  ;; happened to be calculated as happening again a very tiny amount of
  ;; time into the future, so we remove any collisions that involves
  ;; the same two particles (or particle and wall) as last time.
  set collisions filter [item 1 ? != particle1 or
                         item 2 ? != particle2]
  ;; All done.
  set particle1 nobody
  set particle2 nobody

;; check-for-particle-collision is a particle procedure that determines the time it takes
;; to the collision between two particles (if one exists).  It solves for the time by representing
;; the equations of motion for distance, velocity, and time in a quadratic equation of the vector
;; components of the relative velocities and changes in position between the two particles and
;; solves for the time until the next collision

to check-for-particle-collision
  let my-x xcor
  let my-y ycor
  let my-particle-size size
  let my-x-speed speed * dx
  let my-y-speed speed * dy
  ask other particles
    let dpx (xcor - my-x)   ;; relative distance between particles in the x direction
    let dpy (ycor - my-y)    ;; relative distance between particles in the y direction
    let x-speed (speed * dx) ;; speed of other particle in the x direction
    let y-speed (speed * dy) ;; speed of other particle in the x direction
    let dvx (x-speed - my-x-speed) ;; relative speed difference between particles in x direction
    let dvy (y-speed - my-y-speed) ;; relative speed difference between particles in y direction
    let sum-r (((my-particle-size) / 2 ) + (([size] of self) / 2 )) ;; sum of both particle radii

    ;; To figure out what the difference in position (P1) between two particles at a future
    ;; time (t) will be, one would need to know the current difference in position (P0) between the
    ;; two particles and the current difference in the velocity (V0) between the two particles.
    ;; The equation that represents the relationship is:
    ;;   P1 = P0 + t * V0
    ;; we want find when in time (t), P1 would be equal to the sum of both the particle's radii
    ;; (sum-r).  When P1 is equal to is equal to sum-r, the particles will just be touching each
    ;; other at their edges (a single point of contact).
    ;; Therefore we are looking for when:   sum-r =  P0 + t * V0
    ;; This equation is not a simple linear equation, since P0 and V0 should both have x and y
    ;; components in their two dimensional vector representation (calculated as dpx, dpy, and
    ;; dvx, dvy).
    ;; By squaring both sides of the equation, we get:
    ;;   (sum-r) * (sum-r) =  (P0 + t * V0) * (P0 + t * V0)
    ;; When expanded gives:
    ;;   (sum-r ^ 2) = (P0 ^ 2) + (t * PO * V0) + (t * PO * V0) + (t ^ 2 * VO ^ 2)
    ;; Which can be simplified to:
    ;;   0 = (P0 ^ 2) - (sum-r ^ 2) + (2 * PO * V0) * t + (VO ^ 2) * t ^ 2
    ;; Below, we will let p-squared represent:   (P0 ^ 2) - (sum-r ^ 2)
    ;; and pv represent: (2 * PO * V0)
    ;; and v-squared represent: (VO ^ 2)
    ;;  then the equation will simplify to:     0 = p-squared + pv * t + v-squared * t^2

    let p-squared   ((dpx * dpx) + (dpy * dpy)) - (sum-r ^ 2)   ;; p-squared represents difference
    ;; of the square of the radii and the square of the initial positions

    let pv  (2 * ((dpx * dvx) + (dpy * dvy)))  ;; vector product of the position times the velocity
    let v-squared  ((dvx * dvx) + (dvy * dvy)) ;; the square of the difference in speeds
    ;; represented as the sum of the squares of the x-component
    ;; and y-component of relative speeds between the two particles

    ;; p-squared, pv, and v-squared are coefficients in the quadratic equation shown above that
    ;; represents how distance between the particles and relative velocity are related to the time,
    ;; t, at which they will next collide (or when their edges will just be touching)

    ;; Any quadratic equation that is a function of time (t) can be represented as:
    ;;   a*t*t + b*t + c = 0,
    ;; where a, b, and c are the coefficients of the three different terms, and has solutions for t
    ;; that can be found by using the quadratic formula.  The quadratic formula states that if a is
    ;; not 0, then there are two solutions for t, either real or complex.
    ;; t is equal to (b +/- sqrt (b^2 - 4*a*c)) / 2*a
    ;; the portion of this equation that is under a square root is referred to here
    ;; as the determinant, D1.   D1 is equal to (b^2 - 4*a*c)
    ;; and:   a = v-squared, b = pv, and c = p-squared.
    let D1 pv ^ 2 -  (4 * v-squared * p-squared)

    ;; the next test tells us that a collision will happen in the future if
    ;; the determinant, D1 is > 0,  since a positive determinant tells us that there is a
    ;; real solution for the quadratic equation.  Quadratic equations can have solutions
    ;; that are not real (they are square roots of negative numbers).  These are referred
    ;; to as imaginary numbers and for many real world systems that the equations represent
    ;; are not real world states the system can actually end up in.

    ;; Once we determine that a real solution exists, we want to take only one of the two
    ;; possible solutions to the quadratic equation, namely the smaller of the two the solutions:
    ;;  (b - sqrt (b^2 - 4*a*c)) / 2*a
    ;;  which is a solution that represents when the particles first touching on their edges.
    ;;  instead of (b + sqrt (b^2 - 4*a*c)) / 2*a
    ;;  which is a solution that represents a time after the particles have penetrated
    ;;  and are coming back out of each other and when they are just touching on their edges.

    let time-to-collision  -1

    if D1 > 0
      [ set time-to-collision (- pv - sqrt D1) / (2 * v-squared) ]        ;; solution for time step

    ;; if time-to-collision is still -1 there is no collision in the future - no valid solution
    ;; note:  negative values for time-to-collision represent where particles would collide
    ;; if allowed to move backward in time.
    ;; if time-to-collision is greater than 1, then we continue to advance the motion
    ;; of the particles along their current trajectories.  They do not collide yet.

    if time-to-collision > 0
      ;; time-to-collision is relative (ie, a collision will occur one second from now)
      ;; We need to store the absolute time (ie, a collision will occur at time 48.5 seconds.
      ;; So, we add clock to time-to-collision when we store it.
      ;; The entry we add is a three element list of the time to collision and the colliding pair.
      set collisions fput (list (time-to-collision + ticks) self myself)

;; determines when a particle will hit any of the four walls

to check-for-wall-collision  ;; particle procedure
  ;; right & left walls
  let x-speed (speed * dx)
  if x-speed != 0
    [ ;; solve for how long it will take particle to reach right wall
      let right-interval (box-edge - 0.5 - xcor - size / 2) / x-speed
      if right-interval > 0
        [ assign-colliding-wall right-interval "right wall" ]
      ;; solve for time it will take particle to reach left wall
      let left-interval ((- box-edge) + 0.5 - xcor + size / 2) / x-speed
      if left-interval > 0
        [ assign-colliding-wall left-interval "left wall" ] ]
  ;; top & bottom walls
  let y-speed (speed * dy)
  if y-speed != 0
    [ ;; solve for time it will take particle to reach top wall
      let top-interval (box-edge - 0.5 - ycor - size / 2) / y-speed
      if top-interval > 0
        [ assign-colliding-wall top-interval "top wall" ]
      ;; solve for time it will take particle to reach bottom wall
      let bottom-interval ((- box-edge) + 0.5 - ycor + size / 2) / y-speed
      if bottom-interval > 0
        [ assign-colliding-wall bottom-interval "bottom wall" ] ]

to assign-colliding-wall [time-to-collision wall]  ;; particle procedure
  ;; this procedure is used by the check-for-wall-collision procedure
  ;; to assemble the correct particle-wall pair
  ;; time-to-collision is relative (ie, a collision will occur one second from now)
  ;; We need to store the absolute time (ie, a collision will occur at time 48.5 seconds.
  ;; So, we add clock to time-to-collision when we store it.
  let colliding-pair (list (time-to-collision + ticks) self wall)
  set collisions fput colliding-pair collisions

to choose-next-collision
  if collisions = [] [ stop ]
  ;; Sort the list of projected collisions between all the particles into an ordered list.
  ;; Take the smallest time-step from the list (which represents the next collision that will
  ;; happen in time).  Use this time step as the tick-delta for all the particles to move through
  let winner first collisions
  foreach collisions [ if first ? < first winner [ set winner ? ] ]
  ;; winner is now the collision that will occur next
  let dt item 0 winner
  ;; If the next collision is more than 1 in the future,
  ;; only advance the simulation one tick, for smoother animation.
  set tick-delta dt - ticks
  if tick-delta > 1
    [ set tick-delta 1
      set particle1 nobody
      set particle2 nobody
      stop ]
  set particle1 item 1 winner
  set particle2 item 2 winner

to perform-next-collision
  ;; deal with 3 possible cases:
  ;; 1) no collision at all
  if particle1 = nobody [ stop ]
  ;; 2) particle meets wall
  if is-string? particle2
    [ if particle2 = "left wall" or particle2 = "right wall"
        [ ask particle1 [ set heading (- heading) ]
          stop ]
      if particle2 = "top wall" or particle2 = "bottom wall"
        [ ask particle1 [ set heading 180 - heading ]
          stop ] ]
  ;; 3) particle meets particle
  ask particle1 [ collide-with particle2 ]

to collide-with [other-particle]  ;; particle procedure
  ;;; PHASE 1: initial setup
  ;; for convenience, grab some quantities from other-particle
  let mass2 [mass] of other-particle
  let speed2 [speed] of other-particle
  let heading2 [heading] of other-particle
  ;; modified so that theta is heading toward other particle
  let theta towards other-particle

  ;;; 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-particle
  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 particle 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 particle 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 * v1t) + (v1l * v1l))
  ;; 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.
  set energy kinetic-energy

  if v1l != 0 or v1t != 0
    [ set heading (theta - (atan v1l v1t)) ]
  ;; and do the same for other-particle
  ask other-particle [
    set speed sqrt ((v2t ^ 2) + (v2l ^ 2))
    set energy kinetic-energy

    if v2l != 0 or v2t != 0
      [ set heading (theta - (atan v2l v2t)) ]

  ;; PHASE 5: recolor
  ;; since color is based on quantities that may have changed
  ask other-particle [ recolor ]

;;;;;;;;;;;;; particle coloring procedures ;;;;;;;;;;;;;;

to recolor ;; particle procedure
  ;;let avg-speed 1
  ;; avg-speed is assumed to be 0.5, since particles are assigned a random speed between 0 and 1
  ;; particle coloring procedures for visualizing speed with a color palette,
  ;; red are fast particles, blue slow, and green in between.
  ifelse speed < (0.5 * avg-speed) ;; at lower than 50% the average speed
    [ set color blue ]      ;; slow particles colored blue
    [ ifelse speed > (1.5 * avg-speed) ;; above 50% higher the average speed
        [ set color red ]        ;; fast particles colored red
        [ set color green ] ]    ;; medium speed particles colored green

;;;;;;;;;;; time reversal procedure  ;;;;;;;;;;;;;;;;;;;;;

;; Here's a procedure that demonstrates time-reversing the model.
;; You can run it from the command center.  When it finishes,
;; the final particle positions may be slightly different because
;; the amount of time that passes after the reversal might not
;; be exactly the same as the amount that passed before; this
;; doesn't indicate a bug in the model.
;; For larger values of n, you will start to notice larger
;; discrepancies, eventually causing the behavior of the system
;; to diverge totally. Unless the model has some bug we don't know
;; about, this is due to accumulating tiny inaccuracies in the
;; floating point calculations.  Once these inaccuracies accumulate
;; to the point that a collision is missed or an extra collision
;; happens, after that the reversed model will diverge rapidly.

to test-time-reversal [n]
  ask particles [ stamp ]
  while [ticks < n] [ go ]
  let old-ticks ticks
  while [ticks < 2 * old-ticks] [ go ]
  ask particles [ set color white ]

to reverse-time
  ask particles [ rt 180 ]
  set collisions []
  ask particles [ check-for-wall-collision ]
  ask particles [ check-for-particle-collision ]
  ;; the last collision that happened before the model was paused
  ;; (if the model was paused immediately after a collision)
  ;; won't happen again after time is reversed because of the
  ;; "don't do the same collision twice in a row" rule.  We could
  ;; try to fool that rule by setting particle1 and
  ;; particle2 to nobody, but that might not always work,
  ;; because the vagaries of floating point math means that the
  ;; collision might be calculated to be slightly in the past
  ;; (the past that used to be the future!) and be skipped.
  ;; So to be sure, we force the collision to happen:

;;;;;;;;;;;;;;;;;;;;;;;; reporters ;;;;;;;;;;;;;;;;;;;;;;

to-report init-particle-speed
  report 1

to-report max-particle-mass
  report max [mass] of particles

to-report kinetic-energy
   report (0.5 * mass * speed * speed)

to draw-vert-line [ xval ]
  plotxy xval plot-y-min
  plotxy xval plot-y-max

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

There are 15 versions of this model.

Uploaded by When Description Download
Uri Wilensky about 5 years ago Updated to NetLogo 5.0.4 Download this version
Uri Wilensky over 5 years ago Updated version tag Download this version
Uri Wilensky over 5 years ago Updated to version from NetLogo 5.0.3 distribution Download this version
Uri Wilensky over 6 years ago Updated to NetLogo 5.0 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Updated from NetLogo 4.1 Download this version
Uri Wilensky about 8 years ago Model from NetLogo distribution Download this version
Uri Wilensky about 8 years ago Model from NetLogo distribution Download this version
Uri Wilensky about 8 years ago GasLab Circular Particles Download this version

Attached files

File Type Description Last updated
GasLab Circular Particles.png preview Preview for 'GasLab Circular Particles' over 5 years ago, by Uri Wilensky Download

This model does not have any ancestors.

This model does not have any descendants.