Generating an island with stream-power erosion for Godot

Blobber has an island: 121 square km of land that a Python generator builds from one JSON file by uplift and stream-power erosion, with 973 river reaches that all drain to the sea. Godot streams it round the party as 512 m Terrain3D regions. I build the game by directing AI coding agents, and this is what they built on 8 October.

Laying out a knowledge base for world building

Three research agents wrote a knowledge base in parallel before any code, and the generator’s design comes from named pages in it.

It is a folder of Markdown outside the game’s repository, <workspace>/Tools/Blobber/kb_world/: 86 notes and 15 folder READMEs as it stands now, about 121,000 words. Fourteen numbered folders run from 00_policy through world generation, engine limits, ground textures, water and case studies of shipped games to 12_plan. INDEX.md maps a question to the file to start at, and an agent greps the rest.

Every file opens with the same front matter, here the erosion note’s:

source_authority: peer-reviewed papers (high); developers' write-ups and repositories (medium)
software: none (method notes); reference code in C++/OpenGL and Python
version: 1
domain: 01_worldgen/erosion
asset_class: terrain
target_engine: Godot 4.7
license: notes are the project's; code licences as listed
confidence: medium (methods high; costs at 64 km are extrapolations and marked so)

00_policy/README.md ranks sources (papers and official documentation high, blogs medium, forum posts low) and licences: MIT code may be ported into the game, a GPL tool may only be run as a separate program. A *.note.md is the research, each claim tied to a URL and its fetch date, or marked “from memory, not fetched”. A folder’s README.md is the conclusion: options in a table, a recommendation, and numbered work items with their dependencies.

A fourth agent read only the last section of each README and wrote 12_plan/PLAN.md, 547 lines in eight stages, with every number tagged measured, published, estimated or its own arithmetic.

01_worldgen/README.md compares seven ways to make the large shape of land and recommends authored uplift maps with stream-power erosion, in Python and NumPy as a build step beside the game. Its work items W1 to W5 are the next three sections.

Drawing the coast from control points and noise

The Wishbone’s coast is 34 points in tools/worldgen/islands/wishbone.json.

 "coast": [
  [4, 8], [6, 4.5], [10, 2.8], [15, 2.4], [20, 3.2], [24.5, 4.6], [27.8, 7.5], [29.6, 11.5],
  [30, 15], [28.6, 17.6], [28.8, 20.4], [29.6, 23], [28.4, 25.6], [25.8, 26.2], [24, 24.4],
  [24.6, 21.4], [25.6, 18.4], [24.8, 16], [23.6, 13.6], [21.4, 11.4], [18, 10.4], [15, 11],
  [13, 13.2], [12.6, 16.4], [13.8, 19.6], [14.4, 23], [13.6, 26.4], [11.6, 29.4], [8.6, 30.2],
  [6, 28.6], [4.4, 25.4], [3.4, 21.6], [2.4, 17.5], [2.6, 12.5]
 ],

Positions are in units, 32 to a side: a unit is 512 m of the 16,384 m square, one Terrain3D region. Every map is a 1,024 by 1,024 field at 16 m a cell.

A smooth closed curve through the points becomes a signed distance field, positive inland. The field is domain-warped: it is read at positions pushed by two 3-octave fBm fields, then a 5-octave fBm is added, 800 m of amplitude at a 3,400 m wavelength for the Wishbone. A smooth maximum joins on named necks the noise must not cut, and a smooth minimum cuts the straits and a harbor. Land is where the field is over zero.

The noise is gradient noise from an integer hash of lattice points, so a value depends only on seed and position, and two builds are byte for byte identical. Each step is a pure function that writes fields as raw little-endian .bin files with a JSON header, and runs again only when a hash of its inputs changes.

I had six outlines built this way and picked B.

Raising the land and eroding it with the stream power law

The land’s height comes from a map of how fast the ground rises, worn down by the rivers that form on it. Noise alone has no drainage.

The uplift map is drawn from a tectonic story, data in the same file:

 "tectonics": {
  "note": "The inner sea is a rift that the sea has entered. The two arms are its shoulders, tilted up toward the water and falling away gently to the ocean; the Arch is its closed head, older and lower; the Plug's volcano stands on the rift's floor; the Hook is a low block dropped across a fault at the Wrist.",
  "top_m": 1150,
  "base_uplift": 0.065,
  "basement": "schist",
  "lowland": "mudstone",
  "lowland_below": 0.2,
  "valley_floor": "mudstone",
  "ranges": [
   {
    "name": "the Long Arm",
    "why": "the rift's western shoulder: steep to the inner sea, long slopes to the ocean, highest and coldest at its southern end",
    "spine": [[6.2, 11.8], [6.8, 13.4], [7.6, 18], [8.8, 23], [9.8, 27.4]],
    "uplift": 1.0,
    "along": [0.8, 1.0],
    "half_width": 4.4,
    "skew": -0.25,
    "core": "granite",
    "core_width": 0.28
   },

skew steepens one side of a range. Six rock types each carry a hardness and an angle of rest, from granite (2.4, 50 degrees) to mudstone (0.55, 30 degrees).

erosion.py solves the stream power law, dh/dt = u - K * A^0.5 * s: uplift u raises a cell and its river cuts it down by drainage area A and slope s. Each step sends every land cell’s water to the steepest of its eight neighbors, joins each pit to the sea through its lowest pass, and sums area down that tree. Then it solves heights from the sea upward:

def step(lv: Level, h: np.ndarray, k: float, dt: float, m: float, hillslope: float) -> np.ndarray:
    """One implicit step of uplift and stream-power erosion. Returns the new heights."""
    g = lv.grid
    N = g.N
    rt = Routing(g, h)
    area = rt.accumulate(np.append(np.full(N, lv.cell * lv.cell), 0.0))
    pull = (k / lv.hard[:N]) * (area[:N] + hillslope) ** m * dt / rt.length
    water = rt.level[:N] if rt.level is not None else np.full(N, -np.inf)
    base = h[:N] + lv.u[:N] * dt
    reach = lv.tan[:N] * rt.length
    near = rt.length < lv.cell * 1.5
    hn = h.copy()
    recv = rt.recv
    for d in range(1, rt.deepest + 1):
        nodes = rt.nodes_at(d)
        below = np.maximum(hn[recv[nodes]], water[nodes])     # a lake's surface is its bed's base level
        p = pull[nodes]
        b = base[nodes]
        v = (b + p * below) / (1.0 + p)
        v = np.where(near[nodes], np.minimum(v, below + reach[nodes]), v)
        hn[nodes] = v
    lv.hold(hn)
    relax(lv, hn, 1)
    return hn

A cell’s new height uses its receiver’s new height, already solved. That makes the step implicit, so it is stable at a time step of 40. The np.minimum line is talus: no cell stands over its receiver by more than its rock’s angle of rest.

It runs coarse to fine, each grid started from the last:

K is rescaled during the run so the highest point lands on top_m. The valley, the gorge and the passes are held by fields the erosion obeys every step (a height cap, a height floor, kept lakes), and one command checks 18 constraints on the result.

The last pass is Priority-Flood (Barnes 2014). It fills every hollow with 2 mm of fall toward its outlet, then gives back as lakes the hollows of at least 3,000 square m and 1.5 m deep. The generator’s own notes put the erosion at about 11 seconds and the whole island at about a minute.

Streams run straight for hundreds of metres, because each cell has one receiver.

Routing rivers and finding lakes

Rivers, lakes, springs and waterfalls are read off the eroded heights. Nothing is drawn.

A cell is a creek once 0.045 square km drains through it and a river at 0.45. Strahler order sets the width: 1.5, 3, 6, 12, 22 or 40 m. A reach is the channel between two joins, kept as a line of [x, y, bed height] from its head down. This is the first in tables/rivers.json:

  {
   "discharge_m3s": 2.011,
   "drains_km2": 0.1613,
   "drop_m": 2.19,
   "from": "join",
   "id": 0,
   "kind": "creek",
   "length_m": 48.0,
   "meander_m": 0.0,
   "order": 2,
   "points": [
    [
     7384.0,
     1592.0,
     2.19
    ],
    [
     7384.0,
     1576.0,
     1.5
    ],
    [
     7384.0,
     1560.0,
     0.87
    ]
   ],
   "to": "sea",
   "width_m": 3.0
  },

A lake is the cells under the level their hollow spills at. A waterfall candidate is a run of links at least 0.42 steep and 2.2 times their reach’s usual grade that drops 10 m. The tuned land has 90.5 km of river, 208.9 km of creek, 12 lakes and 612 springs.

Was the land flat enough to play on?

The first land wasn’t: 35% of it was under 10 degrees. A second, tuned for play, is 53%.

play.py bands every cell by slope at the angles the game uses: the party slides on ground over 42 degrees. It then runs Dijkstra over a 32 m grid from the starter valley to eight named places, where a step costs its length over the speed the slope allows (9 m/s to 10 degrees, 5 by 25, 2 by 42). A test in Godot drives the party’s real body along those routes with the game’s own keys and counts how often it is slowed, stopped or slides.

The tuned land is a named variant, a block laid over the island file:

    "top_m": 950,
    "base_uplift": 0.04,
    "range_sharpness": 2.3,
    "coast_plain_m": 500,
    "floodplain": {"min_km2": 0.35, "reach_m": 130, "grade": 0.03, "max_fall": 0.09},
    "sea_cliffs": {"shelter_under": 0.62, "inland_m": 110, "share": 0.8, "min_m": 18},
    "valley_ragged": 0.4,
    "valley_roll_m": 12,
    "valley_bar_m": 20,
    "passes": [
     {"name": "the Saddle", "half_width": 1.2, "cut": 0.75},
     {"name": "the Arch Gap", "half_width": 1.2, "cut": 0.75},
     {"name": "the White Gap", "half_width": 1.1, "cut": 0.3}
    ],
    "erosion": {"hillslope_m2": 80000},

range_sharpness (1.6 before) narrows each range and leaves wider foothills. The stream power law only cuts, so its valleys are V-shaped: floodplain lowers the ground beside a river draining 0.35 square km or more to the river’s height plus a 3% rise.

Painting the ground from slope, wetness and aspect

A painter in Python picks the ground textures at every vertex from what the land is doing there.

The first version stamped noise and covered the plain in same-sized brown blobs. Now paint.py works out causes, each a number per vertex at 1 m. Ridge or hollow is height against its mean 48 m round. Wet is exp(-hand / 1.5) times a ramp on the wetness index ln(flow / slope), where hand is height above the nearest drainage. Sunny is how far the slope faces the sun’s bearing.

A biome in biomes.json lists its grounds. The first lies everywhere. Each later one lies where its causes times their amounts, plus noise and a bias, pass a threshold: bare earth on the dry plain weighs ridge 1.2, sunny 0.8 and wet -1.0. Terrain3D blends two textures a vertex, so the painter keeps the strongest two and packs a 32-bit control word:

    base = np.minimum(first, second).astype(np.uint32)
    over = np.maximum(first, second).astype(np.uint32)
    blend = np.where(second > first, share, 1.0 - share)
    control = (base << 27) | (over << 22) | (np.round(blend * 255.0).astype(np.uint32) << 14)

Textures are stored pale and neutral and the terrain’s color map multiplies them, so one grass serves a dry plain and a wet valley. That keeps the count under Terrain3D’s 32 slots. Python can’t write Terrain3D’s region files. It writes raw float32 heights, uint32 control words and RGBA color, and a Godot script wraps each in an Image and saves a Terrain3DRegion. A 512 m region takes 1.16 s to paint and 45 ms to import.

Streaming the terrain round the party

Only the terrain near the party is in memory. The island is 576 regions of 512 m at a vertex a metre.

Terrain3D 1.0.2 loads every region in its data_directory as its node enters the scene tree. TerrainStream empties that property at run time, before the node enters, so the editor still loads and sculpts the whole terrain and the game starts with none. Then, every frame:

	for location in _wanted(point, load_distance):
		if not loaded.has(location) and not _loading.has(location):
			_loading[location] = available[location]
			ResourceLoader.load_threaded_request(available[location], "", false, ResourceLoader.CACHE_MODE_IGNORE)
	for location: Vector2i in _loading.keys():
		var path: String = _loading[location]
		var status := ResourceLoader.load_threaded_get_status(path)
		if status == ResourceLoader.THREAD_LOAD_IN_PROGRESS:
			continue
		_loading.erase(location)
		var region := ResourceLoader.load_threaded_get(path) as Terrain3DRegion
		if region != null and not loaded.has(location):
			_take_in(location, region)
			arrived.append(location)
	var dropped := _drop_far(point)
	if dropped or not arrived.is_empty():
		_tell(arrived)

A region loads on a thread when its nearest edge is within 512 m and is dropped beyond 704 m. _tell calls the plugin’s update_maps once a frame, 2 to 8 ms.

The plugin’s collision is off: its dynamic mode rebuilt everything within 256 m whenever a region came or went, 265 ms each time. The stream builds its own 64 m squares of HeightMapShape3D within 288 m, 1 ms of building a frame. Physics is Jolt, because holes in Terrain3D hang Godot’s built-in physics.

A bake writes one image of the whole world’s heights beside the regions. It answers height queries where nothing is loaded and feeds the map and the far ground: six square grids of 128 by 128 that follow the camera, 16 m between vertices in the first and double in each next, out to 32 km. One view_distance of 12 km sets the fog’s density to 3 over it and the camera’s far plane to 2.5 times it.

Placed things stream in 64 m chunks, built within 160 m and freed beyond 224 m. At the start only what is near is built: 0.3 s for an 8 km test world with 8,000 things. A save keeps only what the party changed, so a new game’s is 51 KB.

Why were the rivers dry at every join?

Every piece of river faded to nothing over 3.5 m at both ends, so two pieces met at a strip of dry bed.

river_stream.gd cuts each reach into pieces of 14 points. A piece becomes a WaterCourse, a Path3D whose script builds a sheet of water fitted to the bed, when it is within 420 m and its ground is loaded. It is freed beyond 640 m.

The agent fixing the water can’t see my screen, so F3 saves the frame and a JSON file beside it:

	"lake": {
		"level": 169.18,
		"metres": 576.2,
		"sheet": 159,
		"text": "lake sheet 159 at level 169.18 m, 576.2 m away"
	},
	"mode": "flying",
	"note": "",
	"number": 1,
	"over_ground": 1.11,
	"pace": "36 m/s",
	"pitch": -14.17,
	"position": [
		-2917.28,
		169.92,
		-2569.49
	],
	"region": [
		-6,
		-6
	],
	"reproduce": "tools/gd.sh shots world/experiments/explore/explore.tscn -- --at=-2917.28,169.92,-2569.49,33.63,-14.17 --hour=9.02 --shot=again.png",
	"river": {
		"has_water": true,
		"kind": "river",
		"metres": 13.5,
		"reach": 382,
		"text": "river reach 382, 12.0 m wide, 13.5 m away, water built",
		"width": 12.0
	},

Each of my 24 pictures came with a reach id and a command that retakes the view with no window. An agent sorted them into seven kinds of fault and fixed the causes:

By the agent’s own count 12 of the 24 are fixed and 9 are better, and in a survey of 64 other spots dry stretches went from 31 to 4. Seams and stepped banks are left.

Why didn’t the river match the lake?

The river’s triangles were wound back to front. Godot turns a back face’s normal round, so the river was lit from underneath and came out darker than the lake under the same shader. The fix is the index order in water_course.gd:

			indices.append_array([first, first + 7, first + 1, first + 1, first + 7, first + 8])

Wound right, the river turned the lake’s pale cyan. I liked the wrong color, so all water is now the blue the river had.

Also done

• One approved model now comes out in five themes by four states of decay: twenty barrels, then twenty orcs on the same 18 bones. The orcs were 59 MB because every texture was stored twice, in the .glb and in a copy Godot wrote out. They are 15 MB.
• The island is three times the size I meant. I’d read Skyrim as 16 square km when it is about 16 square miles (37 to 41 square km), and the Wishbone’s land is 121.
• Seven test areas of the island are painted, one with crimson meat ground.

Next is the island again in a 10.2 km square, about 47 square km of land.