home

author: niplav, created: 2020-11-20, modified: 2024-07-09, language: english, status: in progress, importance: 4, confidence: highly likely

The Diamond-Square algorithm is a terrain-generation algorithm for two dimensions (producing a three-dimensional terrain). I generalize the algorithm to any positive number of dimensions, and analyze the resulting algorithm.

Generalizing the Diamond-Square Algorithm to n Dimensions

Libre de la metáfora y del mito
labra un arduo cristal: el infinito
mapa de Aquel que es todas Sus estrellas.

— Jorge Luis Borges, “Spinoza”, 1964

Introduction

I learned of the diamond-square algorithm by reading through the archives of Chris Wellon's blog, specifically his post on noise fractals and terrain generation. The algorithm is a fairly simple and old one (dating back to the 80s), but not being interested in graphics or game programming, I shelved it as a curiosity.

However, a while later I needed a way to generate high-dimensional landscapes for a simulation, and remembered the algorithm, I felt like I could contribute something here by generalizing the algorithm to produce landscapes in an arbitrary number of dimensions, and that this would be a fun challenge to sharpen my (then fairly weak) Python and numpy skills.

Description

The original (2-dimensional) diamond-square algorithm, in its simplest form, starts with a $2^n+1 \times 2^n+1$ grid of numbers.

It is easiest explained visually:

1. Either a user or the algorithm itself assigns the four corners some values, which can be random.
2. In the diamond step after that, the value in the middle of the grid is determined as the average of the four values in the corners, plus a random value.
3. Next, the middle value every "face" of the grid is determined by the average of the three values in orthogonal directions plus a random value—the square step.
4. The grid is broken down into four sub-grids, and each sub-grid undergoes the diamond step and the square step. The only difference is in the square step: If a point on the grid lies at the face of two sub-grids, it receives the average of all four orthogonal points.
5. The algorithm terminates if each sub-grid is of size $1 \times 1$.

For $n$ dimensions, do that, just higher-dimensional.

We start by initializing an n-dimensional space with zeros, and the corners with random values:

def create_space(dim, size, minval, maxval, factor):
    space=np.zeros([size]*dim)
    corners=(size-1)*get_cornerspos(dim)
    space[*(corners.T)]=np.random.randint(minval, maxval, size=len(corners))

Here, get_cornerspos is just the one-liner return np.array(list(itertools.product([0, 1], repeat=dim))).

We then intialize the variable offsets, and call the recursive diamond-square algorithm:

offsets=[np.array([0]*dim)]
return ndim_diamond_square_rec(space, dim, size, offsets, minval, maxval, factor)

Diamond

Square

Code here. I think this is probably the 2nd-most beautiful code I've ever written, just after this absolute smokeshow.

Analysis

Results

One Dimension

Two Dimensions

Three Dimensions