globals [
  fast average slow         ;; current counts
  avg-speed                 ;; current averages
  clock vsplit vclock       ;; clock variables
  fade-needed?
  aggregate-list  ;; list of the sums of the temperature at each height
]

turtles-own [
  speed mass energy                     ;; turtle info
  v1t v1l tmp-turtle                    ;; collision info (turtle 1)
  heading2 mass2 speed2 v2t v2l turtle2 ;; collision info (turtle 2)
  theta                                 ;; collision info (both turtles)
]

to setup
  ca
  make-box
  set fade-needed? false
  set clock 0
  set vclock 0
  cct number [
    set speed initspeed
    set mass initmass
    if who != 0
      [ random-position ]
    rt random-float 360
    set shape "circle"
    recolor
  ]
  update-variables
  setup-plots
  setup-histograms
  do-plotting
  do-histograms
end

to update-variables
  ask turtles
    [ set energy (0.5 * speed * speed * mass) ]
  set average count turtles with [color = green]
  set slow    count turtles with [color = blue]
  set fast    count turtles with [color = red]
  set avg-speed mean values-from turtles [speed]
  set vsplit ceiling max values-from turtles [speed + gravity]
end

to factor-gravity  ;; turtle procedure
  locals [vx vy]
  set vx (sin heading * speed)
  set vy (cos heading * speed) - gravity / vsplit
  set speed sqrt ((vy * vy) + (vx * vx))
  recolor
  set heading atan vx vy
end

to go
  ask turtles [ bounce ]
  ask turtles [ move ]
  if not any? turtles [ stop ]
  if trace?
    [ ask min-one-of turtles [ who ]
        [ stamp gray
          set fade-needed? true ] ]
  set vclock (vclock + 1)
  if (vclock = vsplit)
  [ set clock (clock + 1)
    set vclock 0
    update-variables
    do-plotting
    do-histograms
    if fade-needed? [ fade-patches ]
  ]
end

to bounce  ;; turtle procedure
  ; if hitting ground or ceiling, reflect heading around y axis
  if pcolor-of patch-ahead 1 = yellow
    [ set heading (180 - heading) ]
end

to move  ;; turtle procedure
  ;; In other GasLab models, we use "jump speed / vsplit" to move the
  ;; turtle the right distance along its current heading.  In this
  ;; model, though, the molecules are affected by gravity as well, so we
  ;; need to offset the turtle vertically by an additional amount.  The
  ;; easiest way to do this is to use "setxy" instead of "jump".
  ;; Trigonometry tells us that "jump speed / vsplit" is equivalent to:
  ;;   setxy (xcor + sin heading * speed / vsplit)
  ;;         (ycor + cos heading * speed / vsplit)
  ;; so to take gravity into account we just need to alter ycor
  ;; by an additional amount given by the classical physics equation:
  ;;   y(t) = 0.5*a*t^2 + v*t + y(t-1)
  ;; but taking vsplit into account, since vsplit is a subdivider of t.
  setxy (xcor + sin heading * speed / vsplit)
        (ycor + cos heading * speed / vsplit - gravity / (2 * vsplit * vsplit))
  factor-gravity
  check-for-collision
end

to check-for-collision  ;; turtle procedure
  if count other-turtles-here = 1
    [ set tmp-turtle random-one-of other-turtles-here
      if ((who > who-of tmp-turtle) and (turtle2 != tmp-turtle))
        [ collide ]
    ]
end

to collide  ;; turtle procedure
  get-turtle2-info
  calculate-velocity-components
  set-new-speed-and-headings
end

to get-turtle2-info  ;; turtle procedure
  set turtle2 tmp-turtle
  set mass2 mass-of turtle2
  set speed2 speed-of turtle2
  set heading2 heading-of turtle2
end

to calculate-velocity-components
  locals [vcm]  ;; CM vel. along dir. theta
  set theta (random-float 360)
  set v1l (speed * sin (theta - heading))
  set v1t (speed * cos (theta - heading))
  set v2l (speed2 * sin (theta - heading2))
  set v2t (speed2 * cos (theta - heading2))
  set vcm (((mass * v1t) + (mass2 * v2t)) / (mass + mass2))
  set v1t (vcm + vcm - v1t)
  set v2t (vcm + vcm - v2t)
end

to set-new-speed-and-headings  ;; turtle procedure
  set speed sqrt ((v1t * v1t) + (v1l * v1l))
  set heading (theta - (atan v1l v1t)) 

  set speed-of turtle2 sqrt ((v2t * v2t) + (v2l * v2l))
  set heading-of turtle2 (theta - (atan v2l v2t))

  recolor
  ask turtle2 [ recolor ]
end

to recolor  ;; turtle procedure
  ifelse speed < (0.5 * initspeed)
    [ set color blue ]
    [ ifelse speed > (1.5 * initspeed)
        [ set color red ]
        [ set color green ] ]
end

to fade-patches
  locals [trace-patches]
  set trace-patches patches with [(pcolor != yellow) and (pcolor != black)]
  ifelse any? trace-patches
    [ ask trace-patches
      [ set pcolor ( pcolor - 0.4 )
        if (not trace?) or (round pcolor = black)
          [ set pcolor black ] ] ]
    [ set fade-needed? false ]
end

to make-box
  ask patches with [abs pycor = (screen-edge-y)]
    [ set pcolor yellow ]
end

to random-position  ;; turtle procedure
  setxy (random-float screen-size-x)
        (screen-edge-y - 1 - random-float (screen-size-y - 2))
end

;;; plotting procedures
to setup-plots
  set-current-plot "Speed Counts"
  set-plot-y-range 0 number

  set-current-plot "Height vs Temp."
  set-plot-x-range 0 (screen-edge-y * 2)
  
  set-current-plot "Aggregate Temp."
  set-plot-x-range 0 (screen-edge-y * 2)
  set aggregate-list fill-list (plot-x-max / 4) 0
end

to do-plotting
  ;; use without-interruption so plots update all at once,
  ;; without intermediate states ever being visible
  without-interruption
    [ set-current-plot "Speed Counts"
      set-current-plot-pen "fast"
      plot fast
      set-current-plot-pen "average"
      plot average
      set-current-plot-pen "slow"
      plot slow
    
      set-current-plot "Height vs Temp."
      plot-pen-reset
      draw-energy-graph
      
      set-current-plot "Aggregate Temp."
      plot-pen-reset
      draw-aggregate-graph aggregate-list ]
end

to setup-histograms
  set-current-plot "Speed Histogram"
  set-plot-x-range 0 (initspeed * 3)
  set-plot-y-range 0 ceiling (number / 10)
  set-current-plot-pen "average"
  set-histogram-num-bars 45
  set-current-plot-pen "slow"
  set-histogram-num-bars 45
  set-current-plot-pen "fast"
  set-histogram-num-bars 45
  set-current-plot-pen "init-avg-speed"
  draw-vert-line initspeed

  set-current-plot "Energy Histogram"
  set-plot-x-range 0 (0.5 * (initspeed * 3) * (initspeed * 3) * initmass)
  set-plot-y-range 0 ceiling (number / 10)
  set-current-plot-pen "average"
  set-histogram-num-bars 45
  set-current-plot-pen "slow"
  set-histogram-num-bars 45
  set-current-plot-pen "fast"
  set-histogram-num-bars 45

  set-current-plot "Density Histogram"
  set-plot-x-range 0 (screen-edge-y * 2)
  set-plot-y-range 0 (screen-edge-x * 2)
  set-histogram-num-bars round (screen-edge-x / 2)
end

to do-histograms
  set-current-plot "Speed Histogram"
  set-current-plot-pen "average"
  histogram-from turtles with [color = green] [speed]
  set-current-plot-pen "slow"
  histogram-from turtles with [color = blue] [speed]
  set-current-plot-pen "fast"
  histogram-from turtles with [color = red] [speed]
  set-current-plot-pen "average-speed"
  plot-pen-reset
  draw-vert-line avg-speed

  set-current-plot "Energy Histogram"
  set-current-plot-pen "average"
  histogram-from turtles with [color = green] [energy]
  set-current-plot-pen "slow"
  histogram-from turtles with [color = blue] [energy]
  set-current-plot-pen "fast"
  histogram-from turtles with [color = red] [energy]

  set-current-plot "Density Histogram"
  histogram-from turtles [ycor + screen-edge-y]
end

to draw-energy-graph
  locals [n x list-pos temperature energy-sum]
  set n (- screen-edge-y)
  repeat (screen-edge-y / 2)
  [ set x turtles with [ (pycor >= n) and (pycor < (n + 4)) ]
    ifelse count x = 0
     [ plot 0 ]
     [ set temperature mean values-from x [ energy ]
       plot temperature
       set list-pos (n + screen-edge-y) / 4
       set aggregate-list (replace-item list-pos aggregate-list ((item list-pos aggregate-list) + temperature))
     ]
    set n n + 4
  ]
end

to draw-aggregate-graph [lst]
  locals [n x]
  set n 0
  repeat length lst
  [ plot item n lst
    set n n + 1
  ]
end

to-report fill-list [n value]
  locals [lst]
  set lst []
  repeat n
  [ set lst (lput value lst) ]
  report lst
end

; draws a vertical line at xval on the current-plot with the current plot-pen
to draw-vert-line [xval]
  plotxy xval plot-y-min
  plot-pen-down
  plotxy xval plot-y-max
  plot-pen-up
end


; *** NetLogo Model Copyright Notice ***
;
; This model was created as part of the project:
; PARTICIPATORY SIMULATIONS: NETWORK-BASED DESIGN FOR SYSTEMS LEARNING IN
; CLASSROOMS.  The project gratefully acknowledges the support of the
; National Science Foundation (REPP program) -- grant number REC #9814682.
;
; Copyright 2002 by Uri Wilensky.  Updated 2002.  All rights reserved.
;
; Permission to use, modify or redistribute this model is hereby granted,
; provided that both of the following requirements are followed:
; a) this copyright notice is included.
; b) this model will not be redistributed for profit without permission
;    from Uri Wilensky.
; Contact Uri Wilensky for appropriate licenses for redistribution for
; profit.
;
; To refer to this model in academic publications, please use:
; Wilensky, U. (2002).  NetLogo GasLab Gravity Box model.
; http://ccl.northwestern.edu/netlogo/models/GasLabGravityBox.
; Center for Connected Learning and Computer-Based Modeling,
; Northwestern University, Evanston, IL.
;
; In other publications, please use:
; Copyright 1998 by Uri Wilensky.  All rights reserved.  See
; http://ccl.northwestern.edu/netlogo/models/GasLabGravityBox
; for terms of use.
;
; *** End of NetLogo Model Copyright Notice ***
@#$#@#$#@
GRAPHICS-WINDOW
285
10
617
363
80
80
2.0
0
10
1
1
1

CC-WINDOW
98
229
276
364
Command Center

MONITOR
216
171
276
220
clock
clock
0
1

MONITOR
9
229
89
278
 avg-speed
avg-speed
2
1

MONITOR
78
171
138
220
average
average
0
1

MONITOR
9
171
69
220
fast
fast
0
1

MONITOR
147
171
207
220
slow
slow
0
1

SWITCH
147
42
276
75
trace?
trace?
0
1
-1000

SLIDER
147
128
276
161
initmass
initmass
1.0
20.0
5.0
1.0
1
NIL

SLIDER
147
85
276
118
number
number
1
2000
400
1
1
NIL

SLIDER
9
128
138
161
initspeed
initspeed
1.0
20.0
5.0
1.0
1
NIL

BUTTON
78
42
138
75
go
go
T
1
T
OBSERVER
NIL

BUTTON
9
42
69
75
NIL
setup
NIL
1
T
OBSERVER
T

PLOT
9
374
216
544
Speed Counts
time
count
0.0
20.0
0.0
400.0
true
false
PENS
"fast" 1.0 0 -65536 true
"average" 1.0 0 -11352576 true
"slow" 1.0 0 -16776961 true

PLOT
225
374
437
544
Speed Histogram
speed
count
0.0
15.0
0.0
40.0
false
false
PENS
"fast" 1.0 1 -65536 true
"average" 1.0 1 -11352576 true
"slow" 1.0 1 -16776961 true
"average-speed" 1.0 0 -7566196 true
"init-avg-speed" 1.0 0 -16777216 true

PLOT
446
374
658
544
Energy Histogram
energy
count
0.0
450.0
0.0
40.0
false
false
PENS
"fast" 1.0 1 -65536 true
"average" 1.0 1 -11352576 true
"slow" 1.0 1 -16776961 true

SLIDER
9
85
138
118
gravity
gravity
0.0
5.0
1.0
0.1
1
NIL

PLOT
667
185
864
341
Height vs Temp.
height
temperature
0.0
160.0
0.0
100.0
true
false
PENS
"default" 4.0 1 -11352576 true

PLOT
667
19
864
175
Density Histogram
height
count
0.0
160.0
0.0
160.0
true
false
PENS
"default" 4.0 1 -11352576 true

PLOT
667
351
864
507
Aggregate Temp.
height
temperature
0.0
160.0
0.0
500.0
true
false
PENS
"default" 4.0 1 -11352576 true

@#$#@#$#@
WHAT IS IT?
-----------
This model is part of the GasLab collection of models.  See GasLab Gas in a Box for an introduction to the GasLab collection.

This model simulates the effect of gravity on gas molecules.  It is very similar to GasLab Atmosphere, but with a "ceiling" over the gas, so molecules can never escape.


HOW IT WORKS
------------
For the basics of how all GasLab models work, see GasLab Gas in a Box.

In this model the molecules behave exactly as in the basic Gas in a Box model except for now gravity acts on the molecules.  Every molecule is given additional velocity downward during each tick, as it would get in a gravitational field.  The molecules bounce off the "ground" as well as the ceiling.


HOW TO USE IT
-------------
Initial settings:
- GRAVITY: strength of the gravitational acceleration
- NUMBER: number of gas molecules
- INITSPEED: initial speed of each molecule
- INITMASS: mass of each molecule

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

Other settings:
- TRACE?: Traces the path of one of the molecules.  This path fades over time to make the screen less cluttered.

Monitors:
- FAST, AVERAGE, SLOW: numbers of molecules with different speeds: fast (red), average (green), and slow (blue).
- LOSSES: number of molecules that have disappeared at the top edge of the screen.
- AVG-SPEED: average speed of the molecules.
- CLOCK: number of ticks that have run.

Plots:
- SPEED COUNTS: plots the number of molecules in each range of speed.
- SPEED HISTOGRAM: speed distribution of all the molecules.  The gray line is the average value, and the black line is the initial average.
- ENERGY HISTOGRAM: distribution of energies of all the molecules, calculated as m*v*v/2.
- HEIGHT VS. TEMPERATURE: shows the temperature of the molecules at each 'layer' of the box.
- DENSITY HISTOGRAM: shows the number of molecules at each 'layer' of the box, i.e. its local density.
- AGGREGATE TEMPERATURE: shows the aggregate sum of the HEIGHT VS. TEMPERATURE plot for the entire model run.


THINGS TO NOTICE
----------------
Try to predict what the screen will look like after a while, and why.

Watch the gray path of one molecule. What can you say about its motion?

Watch the change in density distribution as the model runs.

As the model runs, what happens to the average speed and kinetic energy of the molecules?  If they gain energy, where does it come from?   What happens to the speed and energy distributions?

What is the shape of the path of individual molecules?

What happens to the aggregate temperature plot over time? Is the temperature uniform over the box?


THINGS TO TRY
-------------
What happens when gravity is increased or decreased?

Change the initial number, speed and mass.  What happens to the density distribution?

Does this model come to some sort of equilibrium?  How can you tell when it has been reached?  

Try and find out if the distribution of the molecules in this model is the same as what is predicted by conventional physical laws.  

Try making gravity negative.  

Varying model parameters such as NUMBER and GRAVITY, what do you observe about the aggregate temperature graph?  Can you explain these results?  How would you square these results with a) the fact that "hot air rises" in a room and b) that it is colder as you go higher in elevation?


EXTENDING THE MODEL
-------------------
Try this model with molecules of different masses.  You could color each mass differently to be able to see where they go.  Are their distributions different?  Do the different masses tend to have different average heights?  What does this suggest about the composition of an atmosphere?


This basic model could be used to explore other situations where freely moving molecules have forces on them -- e.g. a centrifuge or charged particles (ions) in an electrical field. 


NETLOGO FEATURES
----------------
Because of the influence of gravity, the molecules follow curved paths.  Since NetLogo models time in discrete steps, these curved paths must be approximated with a series of short straight lines.  This is the source of a slight inaccuracy where the molecules gradually lose energy if the model runs for a long time.  (The effect is as though the collisions with the ground were slightly inelastic.)  The inaccuracy can be reduced by increasing vsplit, but the model will run slower. 


RELATED MODELS
---------------
This model is part of the GasLab suite and curriculum.  See, in particular, Gas in a Box and GasLab Atmosphere.


CREDITS AND REFERENCES
----------------------
To refer to this model in academic publications, please use: Wilensky, U. (1998).  NetLogo GasLab Gravity Box model. http://ccl.northwestern.edu/netlogo/models/GasLabGravityBox. Center for Connected Learning and Computer-Based Modeling, Northwestern University, Evanston, IL.

In other publications, please use: Copyright 1998 by Uri Wilensky.  All rights reserved.  See http://ccl.northwestern.edu/netlogo/models/GasLabGravityBox for terms of use.
@#$#@#$#@
default
true
0
Polygon -7566196 true true 150 5 40 250 150 205 260 250

circle
false
0
Circle -7566196 true true 35 35 230

@#$#@#$#@
NetLogo 2.0beta5
@#$#@#$#@
@#$#@#$#@
@#$#@#$#@