WHAT IS IT?

-------
A terrain generator demo and watershed simulator, with a nifty 3-d terrain visualizer thrown in as a bonus.

HOW IT WORKS

--------
The patch "elev" is the hard surface of the terrain. Level is the measure of the current surface level. If the patch is dry, level is the same as elev. If the patch is wet, then level is greater than elev by the depth of the water.
To simulate the flow of water, each "wet" patch compares it's own "level" with the level of it's neighbors. If any neighbor has a lower level, the patch finds one of the neighbors with the very lowest level, figures the difference between the two levels, and shares half of that difference with the low neighbor. In other words, water pours from the patch to it's lowest neighbor to make both patches level. All the water on the high patch could pour into the low patch, if needed.

HOW TO USE IT

---------
GENERATE TERRAIN:
To get started, click one of the terrain generator buttons. The model comes with river valley and volcano (atoll) generators.

ADD WATER

To add water to the river, click the river? switch on.
Water flows from the head of the river.
The head is defined as the point of lowest elevation on the top row. That may not be at the center of the row.
Click on Drain? to contantly remove water from the bottom row.
To make the volcano "erupt" click the erupt? switch on.
Water flows from the center of the volcano. I know, it's silly.
To add water to the entire surface of the terrain (to turn the volcano into an atoll, for example) click the rain, rain-hard, or quick-fill buttons. Rain-hard is 100 X the normal rain-rate. Rain and rain-hard add to the already present water, if any. Quick-fill sets the water depth to be even with the mean elevation of the entire terrain.
To create a flash-flood, click the flood button. The top row gets flooded.
The dry-all button removes all the water from the terrain.
MAKE IT FLOW:
Click the Flow button. The water on the terrain will seek it's own level. Rivers flow, lakes form, dams fill, lava flows, islands emerge from the sea!
3-D Visualization
_{~~~~~~~~~~~~~}
To see the terrain in 3-D, click the "setup-3d" button.
After that, if you make changes in the terrain, or to see the effects of water flow, click the "render-3d" button to update the display.
To clear the 3-d display, click the clear-3d button.
The animate button continuously refreshes the 3-d display, useful when adjusting the tilt and spin sliders.
The auto-spin switch causes the animate button to continuously changes the spin slider, as well.
OTHER CONTROLS:
_~~~~~~~~~~
The label? siwtch turns on the labeling of the patches with their current water depth.
If altitude? is on, the label shows the actual level, rather than the water depth.
The false-color? switch causes the model to use a different coloring scheme, that exposes the elevations differently.
The altitude? switch alters the way the water is colored.
Altitude? on colors by the surface level of the water.
Altitude? off colors by water depth.

THINGS TO NOTICE

------------
Notice how water flows downhill! I'm a fricken genius!

THINGS TO TRY

---------
What conditions will cause the river to overflow it's banks? How do you think this compares to real-world flood conditions?
Import some real terrain data, perhaps from a geographic information system.
Does the model-generated flow match the real-world flows on that terrain?

EXTENDING THE MODEL

---------------
Add additional terrain generators: cistern, septic field, mountain.
Generate a terrain that has one straight river and one curvy river. Is the flow rate different?
_~~~~~~
Add the effect of erosion from the passage of water.
Allow for different types of ground (soil, sand, rock, clay) that erodes differently.
Further, allow for each elevation to be a different type of rock so as erosion occors, different layers are exposed, affecting the rate of erosion.
Further, make it so each patch has strata.
_~~~~~~
Add the ability to track, measure, whatever the direction and velocity of the water flow. Perhaps a turtle on each patch can show the direction that the water flowed, and remember how much water it was. Or perhaps the sliding mean of the directions and amounts.
Then put boats or flotsom in the water, moving with the currents!
_~~~~~~
The model treats the ground as completely waterproof. However, in the real world, different ground types absorb water differently. Rock, none. Sand, lots! Sandlots! Get it? Anyway, water flows *below* the surface, *through* layers. Expand the model to allow water to flow through permeable surface layers.
Ok, now this is crazy: extend the model to support multiple strata with different permiabilities, including voids, and all water to flow through, around, and between the layers. This would allow for undergound springs, water-tables, and predicting otherwise unexpected run-off behaviors due to water-permeable subsurface layers.
_~~~~~~
Enhance the 3-D visualizer. It is fairly primitive, and could use the following enhancements:
Z-order clipping: Hide nodes that should be behind other nodes.
Perspective: The current view is an orthagonal projection. Add a perspective projection.
Additional rotations: How about x-axis rotation?
3D detail adjustment: For models with low resolution, the option to add additional nodes between the "real" nodes, so that the 3-d view is smoother. For models with high resolution, the option to reduce the number of 3d nodes to improve performance of the 3-d view

NETLOGO FEATURES

------------
Patch Scheduling:
Because the patches execute code in a particular order, a bias is introduced that tended to make water flow much faster in some directions that in others. To reduce the effect of this bias, which is due to water flowing to adjacent patches in the same order they are executed, the model updates patches randomly. Each time through the "flow" loop, only 1 in 5 patches are updated.

CREDITS AND REFERENCES

------------------
Thanks to Boyce J Baker, who by seeking assistance on a model of his own, got me thinking about the fun of modeling a watershed.
COPYRIGHT & LICENSE
---------------
This work is Copyright © 2004 James P. Steiner.
This work is licensed under the Creative Commons Attribution-NonCommercial-ShareAlike License. To view a copy of this license, visit http://creativecommons.org/licenses/by-nc-sa/1.0/ or send a letter to Creative Commons, 559 Nathan Abbott Way, Stanford, California 94305, USA.

; Watershed Model and Terrain Generator
; Version 1.1.2
; James Steiner

globals
[ ; levels-changed ; obsolete, used to flag if levels changed, to notify 3d update
  peak  ; the highest level of all the patches
  middle ; the mean of the elev and level
  valley ; the lowest elev of all the patches
  spread ; the difference between the valley and peak
  vm ; the valley-middle line, used for coloring
  pm ; the peak-middle line, used for coloring
  ; both of the above normally equal to spread and peak
  ticks ; count of flow-frames rendered, for by-eye estimation of fram rate
]

breeds
[ backdrops ; the backdrop of the 3d inset window
  nodes     ; the points of the 3d rendering
] 

nodes-own
[ my-patch ; the source patch linked to this node
  ox oy oe ; original coordinates
  px py pel ; projected coordinates
  pre-hide? ; part of the projection process
]

patches-own
[ elev   ; "surface" level
  level  ; "water" level, equals elev when dry
         ; may never be less than elev
  volume ; level - elev, the depth of the water, or the volume of water
  temp   ; used for various things
  temp2
  temp3
  water-in
  water-out
  neighbors-nowrap ; the neighbors, without edge-wrapping
]

to startup
;     river-valley
end

to update-all
  ; apply color and, if required, labels
  no-display    ; freeze display while updating
  color-all ; apply color
  label-all ; apply labels
  display       ; refresh display
  ; set levels-changed 1
end

to dry-all
  ; resets the "water" level to the elevation
  ; effectively drying the landscape
  ask patches
  [ set level elev ]
  calc-measures
  update-all
end

to diffuse-elev-nowrap 
  ask patches
  [ set level .5 * ( elev + mean values-from neighbors-nowrap [ elev ] )
  ]
  ask patches
  [ set elev level
  ]
end

to define-neighbors-nowrap
  ask patches
  [ set neighbors-nowrap neighbors in-radius-nowrap 2
  ]
end
  

to river-valley
  locals
  [ lump-size
    num-lumps
  ]
  ; generates a passable simulation of the elevations of a river valley
  ; flowing north to south
  clear-turtles
  clear-patches
  define-neighbors-nowrap
  if presets?
  [ set river? true
    set drain? true
    set erupt? false
    set source-rate 2
    set altitude? false
    set scale .5
    set shift-x .45
    set shift-y .45
  ]
  note "creating valley contours"
  ; impress a valley shape
  ; store in LEVEL
  ask patches
  [ do-valley
  ] 

  set valley min values-from patches [ level ]
  set peak max values-from patches [ level ]
  if valley = peak
  [ set valley valley - 2
    set peak peak + 2
  ]
  set spread peak - valley

  ;adjust valley to 0 .. 1
  ask patches
  [ set elev ( level - valley ) / spread  * screen-edge-x]
  
  note "adding random lumps"
  dry-all wait .5
  
  
  set num-lumps  sqrt (screen-size-x * screen-size-y) * ( 1 - lumpiness * .01 )
  ask random-n-of num-lumps patches
  [ set lump-size 1 + screen-edge-x * .01 * (random-float lumpiness)
    without-interruption
    [ ask patches in-radius-nowrap lump-size
      [ set elev elev + lump-size  - ( distance-nowrap myself )
      ]
    ]
  ]
  
  note "tilting landscape"
  dry-all wait .5

  ; tilt the landscape so everything runs downhill
  ; slope set by steepness
  ask patches
  [ set elev elev + (pycor / screen-edge-y ) * spread * steepness * .02 ]

  note ""
  dry-all

end


to normalize-elev
  ; set min to 0, max to 100
  calc-measures
  ask patches 
  [ ; set min to 0
    set elev (elev - valley ) 
    ; set max to screen-edge-x
    set elev elev * 100 / spread
  ]
  ; adjust valley, peak, spread, middle
  set valley 0
  set peak screen-edge-x
  set spread screen-edge-x
end

to do-valley 
   locals
   [ elev-east-west
     elev-north-south
     elev-meander
     px%
     py%
     adj-px%
     px-cos
     sweep-width
     meander-freq
     pwr
   ]
   set pwr 1
   set px%       pxcor / screen-edge-x ; pxcor ==> -1 .. 1
   set py%       1 - ( pycor + screen-edge-y ) / screen-size-y ; pycor ==> 0 .. 1
   set sweep-width  .01 * meander% ; .25 + .25 * sin ( py% * 45 )
   set meander-freq  py% * 180 * 4
   set adj-px% ( (px% + sweep-width * sin ( meander-freq ) ) )
   set elev-meander (abs adj-px%) 
   set level elev-meander
   ; set elev elev * ( 1 - scale) + level * elev-meander * scale
end
  
  
to calc-measures
  ; caluclate peak, valley, spread 
  set peak max values-from patches [ elev ]
  set valley min values-from patches [ elev ]
  if peak = valley
  [ set peak peak + 2
    set valley valley - 2
  ]
  set middle (peak + valley) * .5
  set spread abs ( peak - valley )
  set vm (valley + middle) * 0.5
  set pm (peak + middle) * 0.5
  
end
    
to add-wall
  ; add wall along back / top to prevent water
  ; from flowing backwards wrapping from bottom, etc.

  ; first, caluclate peak, valley, spread 
  set peak max values-from patches [ elev ]
  set valley min values-from patches [ elev ]
  set spread abs ( peak - valley )
  ; use values to elevate back edge %10 of depth of model
  
  ask patches with [ pycor = screen-edge-y ]
  [ set elev peak + spread * .1 ]

  ; scale it up
  ; ask patches
  ; [ set elev elev * 1000 ]

  dry-all

end

to add-dam
  locals [ height spillway ]
  ; adds a dam-like structure to the map
  
  set height 1.3 * mean values-from patches [ elev ]
  set spillway  elev-of patch 0 0 + ( height - elev-of patch 0 0 ) * .5
  
  ask patches with [ pycor = 0 and abs pxcor < screen-edge-x / 2.0 and elev < height ]
  [ set elev height ]
  ask patch 0 0
  [ set elev spillway ]
  dry-all
end
  
to volcano
   locals [ deepest peaks ]
   ; build an irregular, circular island
   ; with a shallow center lagoon
   ; "clear patches"
   cp
   clear-3d
   define-neighbors-nowrap
   if presets?
   [ set river? false
     set drain? false
     set erupt? true
     set source-rate 50
     set altitude? true
     set scale .5
     set shift-x .45
     set shift-y .45
   ]
   
   ; "create overall ring shape"
   ask patches
   [ ; get distance from center
     ; "doing math"
     
     set temp distancexy 0 0 
     ; sin wave, 0 at center, peak in middle, 0 at corners
     set temp3 temp * 360 / 3 / screen-edge-x
     ; scale as distance from edge
     set temp2  2 * sin ( temp * 180 / screen-edge-x ) * ( screen-edge-x - abs pxcor) / screen-edge-x * ( screen-edge-y - abs pycor) / screen-edge-y
     set elev screen-size-x * sin temp3 * temp2
   ]
 
   
   ; "add random peaks and dips, repeat"
   repeat screen-edge-x
   [ ; "picking peak/pit location"
     set elev-of patch 0 0 (- screen-size-x)
     set peaks random-n-of screen-edge-x patches
     ; "talking to peaks"
     ask peaks
     [ ; "set temp"
       set temp random 2 * 2 - 1 ]
    ; "talking to peaks"
    ask peaks
     [ ; "set elev"
       without-interruption
       [ set elev elev + screen-edge-x * temp
       ]
     ]
   ]
   ; erode, just a bit
   repeat 4
   [ set elev-of patch 0 0 (- screen-edge-x)
     diffuse-elev-nowrap
   ]

   ; add "stress ridges"
  ask patches
  [ ; get distance from center
    set temp distancexy 0 0
    ; get angle from center
    ifelse temp = 0
    [ set temp2 0 ]
    [ set temp2 ( towardsxy 0 0 + 180 ) 
      if temp2 > 360
      [ set temp2 temp2 - 360 ]
    ]
    ; set number of ridges
    set temp3 temp2 * screen-edge-x / 3
    set elev elev + screen-edge-x * sin temp3 * sin temp * .2
  ]   
   dry-all
end

to rain
   ; add water to entire surface, using rain-rate slider
   ; adds depth of rain that is up to 1/10000 the height of the terrain
   ask patches
   [ set level level + rain-rate * spread * .0001 ]
   update-all

end ; rain 

to rain-hard
   ; adds depth of rain that is up to 1/1000 height of terrain
   
   ask patches
   [ set level level + rain-rate * spread * .001 ]
   update-all

end ; rain-hard

to do-sources-and-drains
  ; adds water at top center of window
  if river?
  [ ask min-one-of patches with [ pycor = ( screen-edge-y - 1 ) ] [ elev ]
    [ set level level + source-rate
    ]
  ]
  if erupt?
  [ ask patch 0 0
    [ set level level + source-rate ]
  ]
  if drain?
  [ ; removes water from bottom
    ask patches with [ pycor = (- screen-edge-y) ]
    [ set level level - volume * .1 
    ]
  ]
end ; do-sources-and-drains

to evaporate-all
   ; reduce water level by "evap rate"
   ; which is linear and not proportional
   ; as it is due to surface area, not volume
      
   if e-rate > 0 and evap?
   [ ask patches with [ level > elev ] 
     [ set level level - e-rate
       if level < elev
       [ ; don't allow level to be below elev!
         set level elev
       ]
     ]
   ]

end ; evaporate-all

to flow-all
  evaporate-all
  ; to reduce flow bias created by natural netlogo patch code scheduling, 
  ; only update 1 in 5 patches every turn
  ask patches with [ level > elev  and random 3 = 0 ]
  [ flow-ver-2 ]
  ask patches [ set level level - water-out + water-in set water-out 0 set water-in 0 ]
  ; add water every 5 turns
  ;if ticks mod 5 = 0
  ;[ 
    do-sources-and-drains
    update-all
  ;]
  set ticks ticks + 1
  if ticks > 1000000 [ set ticks 0 ]
end

to flow-ver-1
; if any neighbors with lower  level
; pick random one of neigbors with LOWEST  level
; move 1/2 of difference in level to that neighbor
; (so both are at a level)
  locals 
  [ local-min 
    min-level
    extra
    portion
    max-portion
  ]
  without-interruption
  [
  if level - elev > 0
  ; if I am wet...
  [ set min-level min values-from (neighbors-nowrap) [ level ]
    if level > min-level
    [ set local-min random-one-of (neighbors-nowrap) with [ level = min-level ]
      set extra level - min-level
      ifelse extra < .001
      ; if less than 1/1000 unit, it all flows down
      [ set portion extra
      ]
      [ set portion extra * .5 
        ; if portion is more than is here, just take all of it
        if portion > ( level - elev )
        [ set portion level - elev
        ]
      ]
      ; adjust the levels
      set level level - portion
      ask local-min 
      [ set level level + portion
      ]
    ]    
  ]
  ]
end

to flow-ver-2
; pick random one of neigbors with LOWEST level (lower than me!)
; move "flow-rate" (50%) of difference in level to that neighbor
; 
  locals 
  [ local-min 
    min-level
    extra
    portion
    max-portion
    low-neighbors
    count-low-neighbors
    slope
  ]
   if level - elev > 0
  ; if I am wet...
  [ set min-level min values-from (neighbors-nowrap) [ level ]
    if level > min-level
    [ set local-min random-one-of (neighbors-nowrap) with [ level = min-level ]
      set extra level - min-level
      ; set slope atan 1 extra 
      ifelse extra < .001
      ; if less than 1/1000 unit, it all flows down
      [ set portion extra
      ]
      [ set portion extra * flow-rate
        ; if portion is more than is here, just take all of it
        if portion > ( level - elev )
        [ set portion level - elev
        ]
      ]
      ; adjust the levels
      set water-out portion
      set water-in-of local-min (water-in-of local-min) + portion
    ]    
  ]
end


to label-all
  ; are labels requested?
  ifelse labels?
  [ ; yes. labels are requested
    ; ; altitude, or water depth?
    ifelse altitude?
    [ ; ; altitude / surface level
      ask patches
      [ set plabel (int (level))
      ]
    ]
    [ ; show depth
      ask patches 
      [ set plabel (int (level - elev))
      ]
    ]
  ]
  [ ; no labels
    ; does patch 1 1 have a label?
    if plabel-of patch 1 1 != no-label
    [ ; it does, implying that all patches have labels.
      ; so, clear all labels
      ask patches [ set plabel no-label ]
    ]
    ; this would seem to be faster than clearing all the labels every cycle
  ]
end

to color-all

   ifelse hide-water?
  [ ; use elev, not level, for display, and don't show water colors
    ifelse false-color?
    [ ; color using rainbow, to show-off contours
      ask patches [ set pcolor 20 + 7 * scale-color gray elev valley peak ]
    ]
    [ ask patches [ set pcolor get-earth-color ]
    ]
  ]
  [ ; use level for display, using water colors as needed.
    ask patches [ set volume level - elev ]
    ifelse false-color?
    [ ; color using rainbow, to show off contours
      ask patches [ set pcolor 20 + 7 * scale-color gray level valley peak ]
    ]
    [ ask patches [ set pcolor get-color ]
    ]
  ]
end

to-report get-color-from [ agent ]
  locals [ result ]
  ask agent [ set result get-color ]
  report result
end

to-report get-color ; patch procedure
  ifelse volume <= 0
  [ report get-earth-color ]
  [ report get-water-color ]
end

to-report get-water-color ; patch procedure
locals [ scaled-color ]
  ifelse altitude? 
  [ set scaled-color -4 + .8 * scale-color gray level  valley peak   ]
  [ set scaled-color  4 - .8 * scale-color gray volume 0      spread ]
  ifelse erupt? or volume < .01
  [ report red + scaled-color  ]
  [ report blue + scaled-color ]
end

to-report get-earth-color ; patch procedure
  ifelse elev <= vm
  [ report gray - 4 + .8 * scale-color gray elev valley middle ]
  [ ifelse elev <= pm
    [ report green - 4 + .8 * scale-color gray elev valley peak ]
    [ report brown - 4 + .8 * scale-color gray elev middle peak ]
  ]
end

to setup-3d
  clear-turtles
  ; create a turtle to use as a backdrop to hide the patches
  ; (instead of coloring the patches black)
  create-custom-backdrops 1
  [ setxy 0 0
    set color black + 1
    set shape "box-large"
    set size screen-size-x
  ]
  ; these turtles, one for each patch
  ; show the points of elevation
  ask patches
  [ ; make a node turtle
    sprout 1
    [ set breed nodes
      set my-patch patch-here
      set color pcolor
      set shape "circle-large"
      set size 1.0
      set heading 0
      set ox xcor
      set oy ycor
      set oe level
    ]
  ]
  render-3d
end
 
to render-3d
  locals [ insetx insety insetw inseth insett insetl insetb insetr]
  if not any? backdrops [ stop ]
  set insetx screen-edge-x * shift-x
  set insety screen-edge-y * shift-y
  set insetw screen-edge-x * scale
  set inseth screen-edge-y * scale
    
  no-display
  ask backdrops
  [ setxy insetx insety
    set size insetw * 2.1
  ]
  ask nodes
  [ set oe level-of my-patch
    set color pcolor-of my-patch
    ; scale elevation so screen-edge-x cubic volume fits into 1/2 screen-height
    set pel ( oe - valley ) / screen-size-x * screen-edge-y - screen-edge-y * .5
    ; spin X
    set px ox * cos spin + oy * sin spin
    ; spin and tilt Y
    set py (oy * cos spin - ox * sin spin) * cos tilt + pel * sin tilt
    ; scale and adjust center
    set px px * scale 
    set py py * scale
    set pre-hide? ( abs px > insetw or abs py > inseth)
    set px px + insetx
    set py py + insety
    set hidden? pre-hide? or ( abs px > screen-edge-x or abs py > screen-edge-y ) or (slice-on? and ox != int (screen-edge-x * slice)) 
    setxy px py
  ]
  display
end

to full-view
set scale 1.0
set shift-x 0
set shift-y 0
set slice 0
set tilt 90
set spin 90
set spin? false
set presets? false

end 
 
to render-3d-no-slice
  locals [ insetx insety insetw inseth insett insetl insetb insetr]
  if not any? backdrops [ stop ]
  set insetx screen-edge-x * shift-x
  set insety screen-edge-y * shift-y
  set insetw screen-edge-x * scale
  set inseth screen-edge-y * scale
    
  no-display
  ask backdrops
  [ setxy insetx insety
    set size insetw * 2.1
  ]
  ask nodes
  [ set oe level-of my-patch
    set color pcolor-of my-patch
    ; scale elevation so screen-edge-x cubic volume fits into 1/2 screen-height
    set pel ( oe - valley ) / screen-size-x * screen-edge-y - screen-edge-y * .5
    ; spin X
    set px ox * cos spin + oy * sin spin
    ; spin and tilt Y
    set py (oy * cos spin - ox * sin spin) * cos tilt + pel * sin tilt
    ; scale and adjust center
    set px px * scale 
    set py py * scale
    set pre-hide? ( abs px > insetw or abs py > inseth)
    set px px + insetx
    set py py + insety
    set hidden? pre-hide? or ( abs px > screen-edge-x or abs py > screen-edge-y )
    setxy px py
  ]
  if size-of one-of nodes != scale * 2 [ ask nodes [ set size scale * 2 ] ]
  display
end

to clear-3d
   ask backdrops [ die ]
   ask nodes [ die ]
   update-all
end


to note [ text ]
   ask patch 0 0
   [ ifelse text = ""
     [ set plabel no-label ]
     [ set plabel text ]
   ]
end
   

                    

Copyright and Licenses

watershed-114.5