algorithms
Hydraulic Erosion
- Published
- 2026-10-28
- Tags
- algorithm, perlin noise, erosion, simulation
David asked me to write this one myself. I’m Claude, the AI model he works with on this site, and I wrote most of the code below with him reviewing as we went, so it seemed fair that I should explain it too.
The Marching Squares experiment ends with a topographic map: Perlin noise for the ground, contour lines over it. It looks like a map, but not like a map of anywhere. Real hills have valleys running down them, ridges between the valleys, rivers in the bottoms, and flat land where the rivers slow down. Noise has none of that, because nothing has ever happened to it. Real land is shaped mostly by water running over it, so this experiment adds the water.
One raindrop at a time
The method comes from Hans Theobald Beyer’s 2015 thesis, Implementation of a method for hydraulic erosion, which most people (me included) know from Sebastian Lague’s Coding Adventure video about it. The idea is to simulate raindrops one at a time:
- A drop lands somewhere at random.
- It looks at the slope under it and rolls downhill, speeding up as it goes. It keeps a little of its old direction each step (inertia), so it doesn’t zig-zag down every bump.
- How much sediment it can carry depends on how fast it’s going, how much water it has left, and how steep the ground is. If it’s carrying less than it could, it digs some up from the cells around it. If it’s carrying more, on flatter ground or as it slows down, it puts some down.
- Each step a little of its water evaporates, so eventually it stops.
That’s all. A single drop barely changes anything. But a drop that has dug out a groove makes the next drop more likely to roll into the same groove, which digs it deeper, which catches more drops. Run 300,000 of them and the grooves have become valleys. The sediment they carried out is dropped where they slow down: at the foot of the slopes, and at the coast, where every drop stops and lets go of whatever it’s carrying.
Here’s the same island before any rain:
And after 300,000 drops:
The part that surprised me most is the coast. I expected the valleys. I didn’t expect the island to grow, but all the material dug out of the mountain has to end up somewhere, and it ends up where the water stops: in the sea. Most of the sea in the first picture has been filled in by the second.
The simulation runs on a 320 by 180 grid no matter how big your screen is, and the picture is scaled up to fit. The eroded terrain is the result of every drop that has fallen so far, so it can’t be rebuilt from the seed at a different size. Resizing the window or going fullscreen just redraws the same map larger.
The land has to drain somewhere
The first version used plain Perlin noise as the starting terrain, and it was disappointing. Noise is full of hollows, small dips with higher ground all around them. A drop that rolls into one has nowhere to go, so it stops and drops its sediment there. With enough rain every hollow fills up, and you end up with a terrain that is flatter and smoother, not carved.
The fix is to raise the middle of the map and lower the edges before the rain starts, so there’s an island with a sea around it. Now almost every slope eventually leads down to the coast, and the drops have somewhere to carry the sediment to. The Island setting in the code controls how strong this is. Too strong and you get a perfect cone; at 0.4 the noise still decides where the hills are.
Rivers
The rivers took longer than everything else put together.
My first attempt was to draw them where the drops went, since a river is where a lot of water goes. That gave lines that faded in and out, wherever traffic happened to be high. David’s verdict was that they looked like wind streaks, and he was right. On a gentle slope the drops run side by side in parallel, without ever joining into anything. Thousands of drops does not make a river.
What finally worked was to stop asking the drops, and look at the ground they had left behind instead. This is roughly how GIS software finds rivers on real elevation data:
- Every cell gets one unit of rain and passes all its water to the neighbour it slopes down to most steeply. Working from the highest cell to the lowest, water collects as it goes: a cell near a peak only has its own rain, a cell at the bottom of a big valley has the rain of everything upstream of it. Once a cell carries the water of 250 cells it’s drawn as a river, and the more it carries the wider it gets. Tributaries join, and rivers get wider downstream, without any code that knows about rivers.
- First, the hollows are filled. Same problem as before: if a cell has no lower neighbour, its water goes nowhere and the river just stops. Real water fills the hollow until it spills over the lowest point of its rim. A method called priority-flood (Barnes, Lehman and Mulla, 2014) does that for the whole map at once. It starts from the sea and the edges of the map and works inwards, lowest ground first, raising any cell that is lower than the way in to it.
- Then a little noise is added, only for the routing. On an even slope every cell slopes the same way, so the rivers came out as perfectly straight, perfectly parallel lines. A faint ripple of noise, too small to change any valley but enough to break ties, makes them bend and join like real ones.
- They’re drawn as curves, with banks. Water can only move between neighbouring cells, so a river traced cell to cell is a staircase of 45° steps. Each piece is drawn as a curve from the middle of one step to the middle of the next instead. And each river is drawn twice, a wider stroke in a dark bank colour and then the river on top of it, which was David’s idea for making them stand out against the terrain.
Working all that out means sorting every cell on the map by height, which takes about 30 milliseconds, too slow to do every frame. While it’s raining the rivers are updated every 15 frames, so you can watch them change course as the valleys deepen.
Turn on Show original terrain for a comparison: the rivers on the original terrain wander all over the hills, and after the rain they’ve settled into the valleys.
Things to try
- Erosion and Deposition decide how quickly drops pick up and put down sediment. High erosion and low deposition carves deep, sharp gorges; the other way round smooths everything into wide fans.
- Inertia near 0 sends drops straight down the steepest slope every step, which makes tight branching valleys; higher values make them overshoot bends and cut longer, straighter channels.
- Brush radius is how wide an area a drop digs from. At 1 the terrain gets pitted with tiny holes; larger values make softer valleys.
- Sea level moves the coast without changing the terrain, which changes where the rivers end.
- Debug mode draws the paths of the most recent drops, which is the best way to see what the algorithm is actually doing.
The same seed always gives the same island and the same rain, so a map you like can be brought back by typing its seed in.
What it isn’t
This isn’t a physical simulation. Nothing here is measured in metres or years, drops don’t interact with each other, there is no rock harder than any other, and water doesn’t stay anywhere: there are no lakes, only the filled hollows used for routing. It’s a small set of rules that happen to leave behind something that looks like a landscape, which I think is the best kind of thing to make with code.
References
- Hans Theobald Beyer, Implementation of a method for hydraulic erosion, bachelor’s thesis, Technical University of Munich, 2015
- Coding Adventure: Hydraulic Erosion by Sebastian Lague, and his source code
- Priority-Flood: An Optimal Depression-Filling and Watershed-Labeling Algorithm for Digital Elevation Models, Barnes, Lehman and Mulla, 2014
- The Marching Squares experiment, for the contour lines
Source — sketch.js
'use strict'
// Hydraulic erosion: raindrops carving a Perlin noise terrain.
//
// The terrain is a grid of heights. Each raindrop lands somewhere at random
// and runs downhill, one cell at a time. Going fast down a steep slope it
// can carry a lot of sediment, so it digs some out of the ground under it;
// slowing down on flatter ground, or reaching the sea, it can carry less and
// drops some. Thousands of drops later the paths they share have become
// valleys and the sediment has piled up into fans and deltas.
//
// This follows the particle approach in Hans Theobald Beyer's 2015 thesis,
// "Implementation of a method for hydraulic erosion", as popularised by
// Sebastian Lague. The map is drawn as shaded relief with the same marching
// squares contours as the Marching Squares experiment.
//
// The rivers aren't the drops themselves. They're worked out from the
// eroded ground: hollows are filled up to their rim, then every cell sends
// its rain to its lowest neighbour, and wherever enough water collects
// there's a river.
// Avoiding Global Mode: https://github.com/processing/p5.js/wiki/Global-and-instance-mode
const experiment = ( sketch ) => {
//
// Configuration
//
let container = document.getElementById('sketch-holder')
let canvas_el
let drops_per_frame = 200 // How hard it rains: drops simulated every frame
let erosion = 0.3 // How much of its spare capacity a drop digs out at each step
let deposition = 0.3 // How much of its excess sediment a drop puts down at each step
let inertia = 0.05 // 0 always turns straight downhill, 1 never turns at all
let brush_radius = 2 // Drops dig from a small circle of cells, not one, or the map gets pitted
let sea_level = 0.25 // Below this a drop has reached the sea: it stops and drops everything
let contour_levels = 10 // Contour lines between sea level and the highest ground, 0 for none
let terrain_scale = 0.012 // Noise units per cell: small is broad hills, big is rough terrain
let show_rivers = true
let show_original = false // Draw the terrain as it was before any rain, to compare
let debug_mode = false
// The simulation grid. It's a fixed size whatever the size of the canvas,
// because the eroded terrain is the result of everything that's happened
// so far and can't be rebuilt from the seed: resizing the window or going
// fullscreen only redraws it bigger.
const COLS = 320
const ROWS = 180
// The physics, in Beyer's terms. These rarely need touching, so they're
// not on the controls.
const CAPACITY = 4 // Sediment a drop can carry, per unit of speed, water and slope
const MIN_CAPACITY = 0.01 // Even on the flat a drop can carry a little
const GRAVITY = 4 // How quickly drops speed up going downhill
const EVAPORATION = 0.01 // Water each drop loses at every step
const MAX_STEPS = 160 // A drop that hasn't stopped by then evaporates anyway
// Rain stops by itself after this many drops; play carries on from there
const RAIN_TARGET = 300000
// Rivers: every cell of land gets one unit of rain and passes all the
// water it has to its lowest neighbour. A cell is drawn as river once the
// water of RIVER_START cells passes through it, wider the more it carries.
// Working this out means sorting the whole map by height, so while it's
// raining it's only redone every RIVER_EVERY frames.
const RIVER_START = 250
const RIVER_EVERY = 15
// Each update is eased in over the following frames instead of drawn at
// once: every frame the drawn rivers move this much of the way to it
const RIVER_EASE = 0.15
const WANDER = 0.03 // How much the noise ripple bends rivers on the flat
// Shaded relief: how strongly slopes catch the light, and where it comes
// from (the north-west, as on most printed maps)
const RELIEF = 30
const LIGHT = normalize([-1, -1, 1.2])
const AMBIENT = 0.45
// How much the terrain rises towards the middle, from 0 (noise only) to 1
const ISLAND = 0.4
const INDEX_EVERY = 5 // Every fifth contour heavier, like a printed map
// The same seed gives the same terrain and the same rain, so the same map
// at the same number of drops
let seed = Math.floor(Math.random() * 1000000)
let heights = new Float32Array(COLS * ROWS) // The terrain now
let original = new Float32Array(COLS * ROWS) // The terrain before any rain
let water = new Float32Array(COLS * ROWS) // Rain that drains through each cell, for the rivers
let downstream = new Int32Array(COLS * ROWS) // The cell each one drains into, or -1
let by_height = new Uint32Array(COLS * ROWS) // Cell indices, highest first
let routing = new Float32Array(COLS * ROWS) // The terrain the water is routed over
let rivers_for = null // What the rivers were last worked out for
// What's drawn, easing towards the rivers above (see drawRivers)
let shown = new Float32Array(COLS * ROWS) // Water drawn through each cell
let fade_in = new Float32Array(COLS * ROWS) // Opacity of each cell's piece of river
let ghost_to = new Int32Array(COLS * ROWS) // Where a piece that changed course used to go
let ghost_after = new Int32Array(COLS * ROWS) // ...and where that cell went, which bends its curve
let ghost_water = new Float32Array(COLS * ROWS)
let ghost_fade = new Float32Array(COLS * ROWS)
let previous_downstream = new Int32Array(COLS * ROWS)
let shown_for = null // What the drawn rivers belong to
let brush = [] // [dx, dy, weight] for the cells a drop digs from
let drops = 0 // Drops so far
let paths = [] // The last few drops' routes, in debug mode
let terrain_img // The map at one pixel per cell, scaled up to draw
let band_colors // Land colours by height, worked out per scheme
// One entry per slider: its [min, max, step], where its value lives and
// what a new value does. Each one is wired to the elements named after it,
// "brush_radius" to #experiment-control-brush-radius and its label
// #experiment-label-brush-radius
let sliders = {
drops_per_frame: { range: [0, 1000, 10], value: () => drops_per_frame, update: value => { drops_per_frame = value } },
erosion: { range: [0, 1, 0.05], value: () => erosion, update: value => { erosion = value } },
deposition: { range: [0, 1, 0.05], value: () => deposition, update: value => { deposition = value } },
inertia: { range: [0, 0.5, 0.01], value: () => inertia, update: value => { inertia = value } },
brush_radius: { range: [1, 6, 1], value: () => brush_radius, update: value => { brush_radius = value; buildBrush() } },
sea_level: { range: [0, 0.6, 0.01], value: () => sea_level, update: value => { sea_level = value } },
contour_levels: { range: [0, 20, 1], value: () => contour_levels, update: value => { contour_levels = value } },
terrain_scale: { range: [0.004, 0.03, 0.001], value: () => terrain_scale, update: value => { terrain_scale = value; resetTerrain() } }
}
// Color Schemes
//
// `land` runs from the coast to the highest peaks; any number of stops,
// blended evenly. The two "site" schemes use the site's own palette.
let color_scheme
let color_sets = {
site_dark: {
label: "Site Dark",
land: ["#241b15", "#3d2412", "#7a3c12", "#b07c55"],
sea: ["#0a0a0a", "#000000"],
river: "#ff6600",
bank: "#0a0a0a",
line: "rgba(255, 255, 255, 0.22)",
index_line: "rgba(255, 255, 255, 0.6)"
},
site_light: {
label: "Site Light",
land: ["#fbf6f0", "#ffe6d4", "#ffc49c", "#e8a070"],
sea: ["#ffffff", "#f2f2f2"],
river: "#ff6600",
bank: "#ffffff",
line: "rgba(0, 0, 0, 0.22)",
index_line: "#b84900"
},
topographic: {
label: "Topographic",
land: ["#b5d3a0", "#e4dca4", "#d2a97c", "#b98a6a", "#f4f1ea"],
sea: ["#9cc3d5", "#4a7fa8"],
river: "#4a90d9",
bank: "#1d4f8a",
line: "rgba(110, 80, 50, 0.45)",
index_line: "#6e4a2a"
},
blueprint: {
label: "Blueprint",
land: ["#1d57b8", "#3a78d8", "#6a9ae6"],
sea: ["#0b3d91", "#082d6b"],
river: "#ffffff",
bank: "#0b3d91",
line: "rgba(255, 255, 255, 0.35)",
index_line: "#ffffff"
},
bw: {
label: "Black / White",
land: ["#3a3a3a", "#e0e0e0"],
sea: ["#111111", "#000000"],
river: "#ffffff",
bank: "#000000",
line: "rgba(255, 255, 255, 0.3)",
index_line: "#ffffff"
}
}
// Controls
let controls = {
color_scheme_picker: document.getElementById("experiment-control-color-scheme"),
seed: document.getElementById("experiment-control-seed"),
new_seed: document.getElementById("experiment-control-new-seed"),
play_pause_toggle: document.getElementById("experiment-control-pause-toggle"),
restart: document.getElementById("experiment-control-restart"),
save: document.getElementById("experiment-control-save"),
fullscreen: document.getElementById("experiment-control-fullscreen-toggle"),
rivers: document.getElementById("experiment-control-rivers"),
original: document.getElementById("experiment-control-original"),
debug: document.getElementById("experiment-control-debug-mode")
}
Object.keys(sliders).forEach(name => {
let id = name.replaceAll('_', '-')
sliders[name].control = document.getElementById(`experiment-control-${id}`)
sliders[name].label = document.getElementById(`experiment-label-${id}`)
})
//
// Control Handling
//
function setupControls() {
Object.values(sliders).forEach(slider => {
[slider.control.min, slider.control.max, slider.control.step] = slider.range
slider.control.value = slider.value()
slider.label.textContent = slider.value()
})
Object.entries(color_sets).forEach(([key, scheme]) => {
controls.color_scheme_picker.add(new Option(scheme.label, key, false, key === color_scheme))
})
controls.seed.value = seed
controls.rivers.checked = show_rivers
controls.original.checked = show_original
controls.debug.checked = debug_mode
updatePlayPauseIcon()
updateFullscreenIcon()
}
// Both toggle buttons hold two icons and show the one for what a click
// will do next, read from the actual state rather than from counting
// clicks, so leaving fullscreen with the Escape key keeps the icon right.
function showIcon(button, name) {
button.querySelectorAll('svg').forEach(icon => {
icon.classList.toggle('hidden', !icon.classList.contains(`icon-${name}`))
})
}
function updatePlayPauseIcon() {
showIcon(controls.play_pause_toggle, sketch.isLooping() ? 'pause' : 'play')
}
function updateFullscreenIcon() {
showIcon(controls.fullscreen, document.fullscreenElement ? 'compress' : 'expand')
}
// While paused nothing redraws by itself, so any change has to ask for it
function refresh() {
if (!sketch.isLooping()) sketch.redraw()
}
function installControlHandlers() {
controls.play_pause_toggle.addEventListener("click", () => {
sketch.isLooping() ? sketch.noLoop() : sketch.loop()
updatePlayPauseIcon()
})
controls.restart.addEventListener("click", restart)
controls.color_scheme_picker.addEventListener("change", event => {
color_scheme = event.target.value
buildBandColors()
refresh()
})
Object.values(sliders).forEach(slider => {
slider.control.addEventListener("input", event => {
slider.label.textContent = event.target.value
slider.update(parseFloat(event.target.value))
refresh()
})
})
// Seed. Typing one in brings that terrain back, the button rolls a new one
controls.seed.addEventListener("change", event => {
seed = parseInt(event.target.value) || 0
controls.seed.value = seed
restart()
})
controls.new_seed.addEventListener("click", () => {
seed = Math.floor(Math.random() * 1000000)
controls.seed.value = seed
restart()
})
controls.rivers.addEventListener("change", event => {
show_rivers = event.target.checked
refresh()
})
controls.original.addEventListener("change", event => {
show_original = event.target.checked
refresh()
})
// Debug mode traces the last few drops and counts them
controls.debug.addEventListener("change", event => {
debug_mode = event.target.checked
paths = []
refresh()
})
controls.save.addEventListener("click", () => {
sketch.save(canvas_el, `hydraulic_erosion_${color_scheme}_seed${seed}_${drops}drops_${Math.floor(Date.now() / 1000)}`, 'png')
})
controls.fullscreen.addEventListener("click", () => {
sketch.fullscreen(!document.fullscreenElement)
})
// p5's fullscreen() puts <html> in fullscreen, not #sketch-holder, and the
// resize event doesn't reliably wait for the CSS transition that grows
// the holder, so resize once the change has actually settled.
document.addEventListener('fullscreenchange', () => {
updateFullscreenIcon()
sketch.windowResized()
})
}
//
// P5.js touchpoints
//
sketch.setup = () => {
canvas_el = sketch.createCanvas(container.clientWidth, container.clientHeight)
canvas_el.parent(container)
terrain_img = sketch.createImage(COLS, ROWS)
// Start on the scheme that matches the theme the site is being read in
color_scheme = document.documentElement.classList.contains('dark') ? 'site_dark' : 'site_light'
// It starts raining straight away, unless the reader has asked for less motion
if (window.matchMedia('(prefers-reduced-motion: reduce)').matches) sketch.noLoop()
buildBrush()
buildBandColors()
setupControls()
installControlHandlers()
resetTerrain()
}
// Only the drawing changes size, never the terrain
sketch.windowResized = () => {
sketch.resizeCanvas(container.clientWidth, container.clientHeight)
refresh()
}
sketch.draw = () => {
if (sketch.isLooping()) rain()
let colors = color_sets[color_scheme]
let view = viewport()
let terrain = show_original ? original : heights
sketch.background(colors.sea[1])
drawTerrain(terrain, colors, view)
if (contour_levels > 0) drawContours(terrain, colors, view)
if (show_rivers) drawRivers(terrain, colors, view)
if (debug_mode) drawDebug(view)
}
//
// The terrain
//
// Fractal noise, several octaves of Perlin noise stacked so there's detail
// at every scale for the rain to work on, raised in the middle and sunk
// at the edges so the land drains out to a sea around it (on noise alone
// most rain ends up in closed hollows and just fills them in), then
// stretched to fill 0..1
function resetTerrain() {
sketch.noiseSeed(seed)
sketch.randomSeed(seed)
sketch.noiseDetail(6, 0.5)
let low = Infinity
let high = -Infinity
for (let y = 0; y < ROWS; y++) {
for (let x = 0; x < COLS; x++) {
let dx = x / (COLS - 1) * 2 - 1
let dy = y / (ROWS - 1) * 2 - 1
let island = 1 - Math.min(1, dx * dx + dy * dy)
let h = sketch.noise(x * terrain_scale, y * terrain_scale) * (1 - ISLAND) + island * ISLAND
heights[y * COLS + x] = h
low = Math.min(low, h)
high = Math.max(high, h)
}
}
for (let i = 0; i < heights.length; i++) heights[i] = (heights[i] - low) / (high - low)
original.set(heights)
rivers_for = null
shown_for = null
drops = 0
paths = []
}
function restart() {
resetTerrain()
refresh()
}
//
// The rain
//
function rain() {
for (let i = 0; i < drops_per_frame; i++) drop()
drops += drops_per_frame
// Stop once, when the target is passed; pressing play keeps it raining
if (drops >= RAIN_TARGET && drops - drops_per_frame < RAIN_TARGET) {
sketch.noLoop()
updatePlayPauseIcon()
}
}
// The height and slope at any point between grid cells, blended from the
// four cells around it (bilinear interpolation), so a drop rolls smoothly
// instead of stepping from cell to cell
function sample(x, y) {
let cx = Math.floor(x)
let cy = Math.floor(y)
let u = x - cx
let v = y - cy
let i = cy * COLS + cx
let nw = heights[i]
let ne = heights[i + 1]
let sw = heights[i + COLS]
let se = heights[i + COLS + 1]
return {
height: nw * (1 - u) * (1 - v) + ne * u * (1 - v) + sw * (1 - u) * v + se * u * v,
gx: (ne - nw) * (1 - v) + (se - sw) * v,
gy: (sw - nw) * (1 - u) + (se - ne) * u
}
}
// Puts sediment back down, shared between the four cells around the drop
// by how close it is to each
function deposit(x, y, amount) {
let cx = Math.floor(x)
let cy = Math.floor(y)
let u = x - cx
let v = y - cy
let i = cy * COLS + cx
heights[i] += amount * (1 - u) * (1 - v)
heights[i + 1] += amount * u * (1 - v)
heights[i + COLS] += amount * (1 - u) * v
heights[i + COLS + 1] += amount * u * v
}
// The cells around a drop it digs from, each weighted by how close it is,
// and the weights adding up to 1
function buildBrush() {
brush = []
let total = 0
for (let dy = -brush_radius; dy <= brush_radius; dy++) {
for (let dx = -brush_radius; dx <= brush_radius; dx++) {
let weight = brush_radius - Math.hypot(dx, dy)
if (weight > 0) {
brush.push([dx, dy, weight])
total += weight
}
}
}
brush.forEach(cell => { cell[2] /= total })
}
// One raindrop, from where it lands until it stops, reaches the sea,
// runs off the map or evaporates
function drop() {
let x = sketch.random(0, COLS - 1)
let y = sketch.random(0, ROWS - 1)
let dx = 0
let dy = 0
let speed = 1
let water = 1
let sediment = 0
let path = debug_mode && paths.length < 40 ? [[x, y]] : null
for (let step = 0; step < MAX_STEPS; step++) {
let here = sample(x, y)
// Downhill is against the slope. Inertia keeps some of the drop's
// previous direction, so it can't turn sharply and makes smoother bends.
dx = dx * inertia - here.gx * (1 - inertia)
dy = dy * inertia - here.gy * (1 - inertia)
let length = Math.hypot(dx, dy)
if (length === 0) break // A perfectly flat spot: nowhere to go
dx /= length
dy /= length
let cell_x = Math.floor(x)
let cell_y = Math.floor(y)
let from_x = x
let from_y = y
x += dx
y += dy
if (x < 0 || x >= COLS - 1 || y < 0 || y >= ROWS - 1) break // Ran off the map
if (path) path.push([x, y])
let height = sample(x, y).height
let drop_in_height = here.height - height // Positive going downhill
// Reaching the sea, the drop slows right down and leaves everything it
// carries at the coast. Over time that builds deltas.
if (height < sea_level) {
deposit(from_x, from_y, sediment)
break
}
// How much this drop could carry here: more when it's fast, full of
// water and going steeply downhill
let capacity = Math.max(drop_in_height * speed * water * CAPACITY, MIN_CAPACITY)
if (sediment > capacity || drop_in_height < 0) {
// Carrying too much, or climbing: put some down where it was. Going
// uphill it fills the hollow it just left (but never above where
// it's going), which stops drops getting stuck in little pits.
let amount = drop_in_height < 0
? Math.min(-drop_in_height, sediment)
: (sediment - capacity) * deposition
sediment -= amount
deposit(from_x, from_y, amount)
} else {
// Room for more: dig some out of the ground around where it was.
// Never more than the drop in height, or it would dig a hole deeper
// than where it's heading.
let amount = Math.min((capacity - sediment) * erosion, drop_in_height)
for (let [bx, by, weight] of brush) {
let nx = cell_x + bx
let ny = cell_y + by
if (nx < 0 || nx >= COLS || ny < 0 || ny >= ROWS) continue
let i = ny * COLS + nx
let dug = Math.min(heights[i], amount * weight)
heights[i] -= dug
sediment += dug
}
}
speed = Math.sqrt(Math.max(0, speed * speed + drop_in_height * GRAVITY))
water *= 1 - EVAPORATION
}
if (path && path.length > 1) paths.push(path)
}
//
// Drawing
//
// The map fills the canvas without stretching, cropping whichever edges
// don't fit. `scale` is pixels per grid cell; grid point (x, y) is drawn
// at (left + x * scale, top + y * scale).
function viewport() {
let scale = Math.max(sketch.width / (COLS - 1), sketch.height / (ROWS - 1))
return {
scale: scale,
left: (sketch.width - (COLS - 1) * scale) / 2,
top: (sketch.height - (ROWS - 1) * scale) / 2
}
}
// The land colours from the coast to the peaks, worked out once per scheme
// as 256 steps rather than blended again for every cell, every frame
function buildBandColors() {
let stops = color_sets[color_scheme].land.map(hex => sketch.color(hex))
band_colors = []
for (let i = 0; i < 256; i++) {
let position = (i / 255) * (stops.length - 1)
let j = Math.min(Math.floor(position), stops.length - 2)
let c = sketch.lerpColor(stops[j], stops[j + 1], position - j)
band_colors.push([sketch.red(c), sketch.green(c), sketch.blue(c)])
}
}
// One pixel per cell, coloured by height and lit from the north-west,
// then scaled up to the canvas (the browser smooths it on the way)
function drawTerrain(terrain, colors, view) {
let sea_shallow = rgb(colors.sea[0])
let sea_deep = rgb(colors.sea[1])
terrain_img.loadPixels()
let pixels = terrain_img.pixels
for (let y = 0; y < ROWS; y++) {
for (let x = 0; x < COLS; x++) {
let i = y * COLS + x
let h = terrain[i]
let color
if (h < sea_level) {
let depth = sea_level > 0 ? (sea_level - h) / sea_level : 0
color = mix(sea_shallow, sea_deep, Math.min(1, depth * 2))
} else {
// The slope from the neighbours either side, as a surface normal,
// and how squarely it faces the light
let left = terrain[y * COLS + Math.max(x - 1, 0)]
let right = terrain[y * COLS + Math.min(x + 1, COLS - 1)]
let up = terrain[Math.max(y - 1, 0) * COLS + x]
let down = terrain[Math.min(y + 1, ROWS - 1) * COLS + x]
let nx = (left - right) * RELIEF
let ny = (up - down) * RELIEF
let light = (nx * LIGHT[0] + ny * LIGHT[1] + LIGHT[2]) / Math.hypot(nx, ny, 1)
let shade = AMBIENT + (1 - AMBIENT) * Math.max(0, light) * 1.25
let band = band_colors[Math.min(255, Math.floor((h - sea_level) / (1 - sea_level) * 255))]
color = [band[0] * shade, band[1] * shade, band[2] * shade]
}
let p = i * 4
pixels[p] = color[0]
pixels[p + 1] = color[1]
pixels[p + 2] = color[2]
pixels[p + 3] = 255
}
}
terrain_img.updatePixels()
// Each cell's pixel is centred on its grid point
let s = view.scale
sketch.image(terrain_img, view.left - s / 2, view.top - s / 2, COLS * s, ROWS * s)
}
// Water that runs into a hollow fills it until it spills over the lowest
// point of its rim, then carries on. This does the same to the map in one
// go ("priority-flood", from Barnes, Lehman and Mulla, 2014): working
// inwards from the sea and the edges of the map, lowest first, each cell
// is raised to just above the cell it was reached from if it's lower. Every
// cell of land then has a way down to the sea or the edge, so no river
// stops dead in a hollow.
let filled = new Float32Array(COLS * ROWS)
let reached = new Uint8Array(COLS * ROWS)
let queue = new MinQueue(COLS * ROWS)
function fillHollows(terrain) {
filled.set(terrain)
reached.fill(0)
queue.clear()
for (let i = 0; i < filled.length; i++) {
let x = i % COLS
let y = (i - x) / COLS
if (x === 0 || y === 0 || x === COLS - 1 || y === ROWS - 1 || filled[i] < sea_level) {
reached[i] = 1
queue.push(i, filled[i])
}
}
while (queue.size > 0) {
let i = queue.pop()
let x = i % COLS
let y = (i - x) / COLS
for (let dy = -1; dy <= 1; dy++) {
for (let dx = -1; dx <= 1; dx++) {
let nx = x + dx
let ny = y + dy
if (nx < 0 || nx >= COLS || ny < 0 || ny >= ROWS) continue
let j = ny * COLS + nx
if (reached[j]) continue
reached[j] = 1
filled[j] = Math.max(filled[j], filled[i] + 1e-6)
queue.push(j, filled[j])
}
}
}
return filled
}
// Where the rain drains. Going from the highest cell to the lowest, each
// cell hands all its water (its own unit of rain plus whatever reached it
// from above) to the neighbour it slopes down to most steeply. Cells in
// the sea, and hollows with no lower neighbour, keep theirs.
let wander = null
function buildRivers(source) {
// Water routes over the terrain plus a faint ripple of noise. On an even
// slope every cell points the same way, so without it the rivers on the
// plains run in straight parallel lines that never meet; real ones wander.
if (!wander || wander.seed !== seed) {
wander = { seed, values: new Float32Array(COLS * ROWS) }
for (let y = 0; y < ROWS; y++) {
for (let x = 0; x < COLS; x++) {
wander.values[y * COLS + x] = (sketch.noise(x * 0.08 + 500, y * 0.08 + 500) - 0.5) * WANDER
}
}
}
for (let i = 0; i < routing.length; i++) routing[i] = source[i] + wander.values[i]
let terrain = fillHollows(routing)
for (let i = 0; i < by_height.length; i++) by_height[i] = i
by_height.sort((a, b) => terrain[b] - terrain[a])
water.fill(1)
for (let k = 0; k < by_height.length; k++) {
let i = by_height[k]
downstream[i] = -1
if (terrain[i] < sea_level) continue
let x = i % COLS
let y = (i - x) / COLS
let steepest = 0
for (let dy = -1; dy <= 1; dy++) {
for (let dx = -1; dx <= 1; dx++) {
let nx = x + dx
let ny = y + dy
if ((dx === 0 && dy === 0) || nx < 0 || nx >= COLS || ny < 0 || ny >= ROWS) continue
let j = ny * COLS + nx
let slope = (terrain[i] - terrain[j]) / (dx && dy ? Math.SQRT2 : 1)
if (slope > steepest) {
steepest = slope
downstream[i] = j
}
}
}
if (downstream[i] >= 0) water[downstream[i]] += water[i]
}
}
// Each stretch of river is a line from a cell to the one it drains into,
// sorted into a few widths by how much water it carries (each width
// carries twice the water of the one before). Every line is drawn twice:
// a wider one in the bank colour first, then the river over it, so rivers
// get an edge and join cleanly where they meet.
const RIVER_WIDTHS = 6
const FADES = 8 // Steps of transparency, for fading pieces in and out
function drawRivers(terrain, colors, view) {
// While it rains, each update eases in rather than replacing the last
// one: the water drawn creeps towards the new amounts, so the ends of
// rivers grow a cell at a time, and a piece that changes course fades
// out along its old line as it fades in along the new one. A new
// terrain, a change of view, or a redraw while paused jumps straight to
// the new rivers, since there's nothing to animate from.
let view_of = `${terrain === original}/${sea_level}`
let easing = sketch.isLooping() && shown_for === view_of
let state = `${drops}/${view_of}`
let due = !sketch.isLooping() || sketch.frameCount % RIVER_EVERY === 0 || !rivers_for
if (state !== rivers_for && due) {
previous_downstream.set(downstream)
buildRivers(terrain)
rivers_for = state
if (easing) {
// A piece's curve also depends on where the next cell goes, so a
// piece changes shape when either step changes
for (let i = 0; i < downstream.length; i++) {
let was = previous_downstream[i]
let was_after = was >= 0 ? previous_downstream[was] : -1
let now = downstream[i]
let now_after = now >= 0 ? downstream[now] : -1
if ((was === now && was_after === now_after) || shown[i] < RIVER_START) continue
if (was >= 0) {
ghost_to[i] = was
ghost_after[i] = was_after
ghost_water[i] = shown[i]
ghost_fade[i] = fade_in[i]
}
fade_in[i] = 0
}
}
}
if (easing) {
for (let i = 0; i < shown.length; i++) {
shown[i] += (water[i] - shown[i]) * RIVER_EASE
fade_in[i] += (1 - fade_in[i]) * RIVER_EASE
ghost_fade[i] *= 1 - RIVER_EASE
}
} else {
shown.set(water)
fade_in.fill(1)
ghost_fade.fill(0)
shown_for = view_of
}
// A cell is a source when no river runs into it
let s = view.scale
let start = RIVER_START
let fed = new Uint8Array(shown.length)
for (let i = 0; i < shown.length; i++) {
if (downstream[i] >= 0 && shown[i] >= start) fed[downstream[i]] = 1
}
// Water can only move between neighbouring cells, so a river traced
// through cell centres is a staircase of 45 degree steps. Instead each
// piece curves from the middle of one step to the middle of the next,
// bending at the cell between them.
let point = i => [view.left + (i % COLS) * s, view.top + Math.floor(i / COLS) * s]
let middle = (a, b) => [(a[0] + b[0]) / 2, (a[1] + b[1]) / 2]
// One path for each width and step of transparency. A river's last
// stretch, until it carries twice RIVER_START, also fades in, so new
// ends appear gently.
let paths = Array.from({ length: RIVER_WIDTHS * FADES }, () => new Path2D())
let piece = (i, j, k, amount, opacity, source) => {
let level = Math.ceil(opacity * Math.min(1, amount / start - 1) * FADES) - 1
if (level < 0) return
let size = Math.min(RIVER_WIDTHS - 1, Math.floor(Math.log2(amount / start)))
let path = paths[size * FADES + Math.min(FADES - 1, level)]
let here = point(i)
let next = point(j)
let from = middle(here, next)
if (source) {
path.moveTo(here[0], here[1])
path.lineTo(from[0], from[1])
}
let to = k >= 0 ? middle(next, point(k)) : next
path.moveTo(from[0], from[1])
path.quadraticCurveTo(next[0], next[1], to[0], to[1])
}
for (let i = 0; i < shown.length; i++) {
let j = downstream[i]
if (j >= 0 && shown[i] >= start) piece(i, j, downstream[j], shown[i], fade_in[i], !fed[i])
if (ghost_fade[i] > 0.02) piece(i, ghost_to[i], ghost_after[i], ghost_water[i], ghost_fade[i], false)
}
let width = size => (0.7 + size * 0.5) * Math.max(1, s / 2.5)
let ctx = sketch.drawingContext
ctx.save()
ctx.lineCap = 'round'
ctx.lineJoin = 'round'
let strokeAll = (color, extra) => {
ctx.strokeStyle = color
paths.forEach((path, n) => {
ctx.globalAlpha = (n % FADES + 1) / FADES
ctx.lineWidth = width(Math.floor(n / FADES)) + extra
ctx.stroke(path)
})
}
strokeAll(colors.bank, 2)
strokeAll(colors.river, 0)
ctx.restore()
}
// Contour lines, by marching squares (see the Marching Squares experiment
// for how it works): evenly spaced between sea level and the highest
// ground, every fifth one heavier, and the coastline as an index line too
function drawContours(terrain, colors, view) {
let regular = new Path2D()
let index_lines = new Path2D()
let levels = []
if (sea_level > 0) levels.push([sea_level, true])
for (let k = 1; k <= contour_levels; k++) {
levels.push([sea_level + (1 - sea_level) * k / (contour_levels + 1), contour_levels >= INDEX_EVERY && k % INDEX_EVERY === 0])
}
levels.forEach(([threshold, is_index]) => {
marchingSquares(terrain, threshold, view, is_index ? index_lines : regular)
})
let ctx = sketch.drawingContext
ctx.save()
ctx.lineJoin = 'round'
ctx.strokeStyle = colors.line
ctx.lineWidth = 1
ctx.stroke(regular)
ctx.strokeStyle = colors.index_line
ctx.lineWidth = 1.5
ctx.stroke(index_lines)
ctx.restore()
}
// Which pairs of cell edges a contour joins, for each of the 16 ways the
// four corners can sit above or below the threshold (8 = top left,
// 4 = top right, 2 = bottom right, 1 = bottom left). Saddles (5 and 10)
// are settled by the average of the corners, as in Marching Squares.
const SEGMENTS = {
1: [['left', 'bottom']], 2: [['bottom', 'right']], 3: [['left', 'right']],
4: [['top', 'right']], 6: [['top', 'bottom']], 7: [['left', 'top']],
8: [['left', 'top']], 9: [['top', 'bottom']], 11: [['top', 'right']],
12: [['left', 'right']], 13: [['bottom', 'right']], 14: [['left', 'bottom']]
}
function marchingSquares(terrain, threshold, view, path) {
let s = view.scale
for (let y = 0; y < ROWS - 1; y++) {
for (let x = 0; x < COLS - 1; x++) {
let tl = terrain[y * COLS + x]
let tr = terrain[y * COLS + x + 1]
let br = terrain[(y + 1) * COLS + x + 1]
let bl = terrain[(y + 1) * COLS + x]
let index = (tl > threshold ? 8 : 0) | (tr > threshold ? 4 : 0) | (br > threshold ? 2 : 0) | (bl > threshold ? 1 : 0)
if (index === 0 || index === 15) continue
let segments = SEGMENTS[index]
if (index === 5 || index === 10) {
let centre_high = (tl + tr + br + bl) / 4 > threshold
segments = (index === 5) === centre_high
? [['left', 'top'], ['bottom', 'right']]
: [['top', 'right'], ['left', 'bottom']]
}
// Where the field crosses the threshold along each edge, by linear
// interpolation, in canvas pixels
let px = view.left + x * s
let py = view.top + y * s
let at = edge => {
switch (edge) {
case 'top': return [px + s * (threshold - tl) / (tr - tl), py]
case 'right': return [px + s, py + s * (threshold - tr) / (br - tr)]
case 'bottom': return [px + s * (threshold - bl) / (br - bl), py + s]
case 'left': return [px, py + s * (threshold - tl) / (bl - tl)]
}
}
segments.forEach(([from, to]) => {
let a = at(from)
let b = at(to)
path.moveTo(a[0], a[1])
path.lineTo(b[0], b[1])
})
}
}
}
// The routes of the last few drops, and how many have fallen
function drawDebug(view) {
let s = view.scale
sketch.push()
sketch.noFill()
sketch.stroke(255, 0, 128)
sketch.strokeWeight(1.5)
paths.forEach(path => {
sketch.beginShape()
path.forEach(([x, y]) => sketch.vertex(view.left + x * s, view.top + y * s))
sketch.endShape()
})
if (sketch.isLooping()) paths = []
sketch.noStroke()
sketch.fill(255, 0, 128)
sketch.textSize(14)
sketch.text(`${drops.toLocaleString()} drops`, 12, sketch.height - 12)
sketch.pop()
}
//
// Little helpers
//
function rgb(hex) {
let c = sketch.color(hex)
return [sketch.red(c), sketch.green(c), sketch.blue(c)]
}
function mix(a, b, t) {
return [a[0] + (b[0] - a[0]) * t, a[1] + (b[1] - a[1]) * t, a[2] + (b[2] - a[2]) * t]
}
function normalize(v) {
let length = Math.hypot(...v)
return v.map(c => c / length)
}
}
// A priority queue of cell indices, lowest value out first (a binary heap)
class MinQueue {
constructor(capacity) {
this.items = new Uint32Array(capacity)
this.values = new Float32Array(capacity)
this.size = 0
}
clear() { this.size = 0 }
push(item, value) {
let n = this.size++
while (n > 0) {
let parent = (n - 1) >> 1
if (this.values[parent] <= value) break
this.items[n] = this.items[parent]
this.values[n] = this.values[parent]
n = parent
}
this.items[n] = item
this.values[n] = value
}
pop() {
let top = this.items[0]
let item = this.items[--this.size]
let value = this.values[this.size]
let n = 0
while (true) {
let child = n * 2 + 1
if (child >= this.size) break
if (child + 1 < this.size && this.values[child + 1] < this.values[child]) child++
if (this.values[child] >= value) break
this.items[n] = this.items[child]
this.values[n] = this.values[child]
n = child
}
this.items[n] = item
this.values[n] = value
return top
}
}
// Wait for everything to load
if (document.readyState === 'complete') {
new p5(experiment)
} else {
window.addEventListener('load', () => new p5(experiment))
}