Saturday, October 27, 2012

Creating an Earth

A while ago I decided I wanted to create something that looks like the surface of a planet, complete with continents & oceans and all. Since I've only been on a small handful of planets, I decided that I'd approximate this by creating something like the Earth on the computer (without cheating and just copying the real Earth). 

Where should I start? Well, let's see what the facts we know about the Earth tell us about how to create a new planet on the computer.

Observation 1: Looking at a map of the Earth, we only see the heights of the surface.

So let's describe just the heights of the Earth's surface.

Observation 2 : The Earth is a sphere.

So (wait for it) we need to describe the height on a spherical surface.

Now we can recast our problem of making an Earth more precisely mathematically. We want to know the heights of the planet's surface at each point on the Earth. So we're looking for field (the height of the planet) defined on the surface of a sphere (the different spots on the planet).

Just like a function on the real line can be expanded in terms of its Fourier components, almost any function on the surface of a sphere can be expanded as a sum of spherical harmonics Ylm. This means we can write the height h of our planets surfaces as

h(\theta, \phi) = \sum A_{lm}Y_l^m(\theta, \phi) \quad (1)
If we figure out what the coefficients A of the sum should be, then we can start making some Earths! Let's see if we can use some other facts about the Earth's surface to get get a handle on what coefficients to use.

Observation 3: I don't know every detail of the Earth's surface, whose history is impossibly complicated.

I'll capture this lack-of-knowledge by describing the surface of our imaginary planet as some sort of random variable. Equation (1) suggests that we can do this by making the coefficients A random variables. At some point we need to make an executive decision on what type of random variable we'll use. For various reasons, [1] I decided I'd use a Gaussian random variable with mean 0 and standard deviation alm:
A_{lm} = a_{lm} N(0,1)
(Here I'm using the notation that N(m,v) is a normal or Gaussian random variable with mean m and variance v. If you multiply a Gaussian random variable by a constant a, it's the same as multiplying the variance by a^2, so a*N(0,1) and N(0,a^2) are the same thing.)

Observation 4: The heights of the surface of the Earth are more-or-less independent of their position on the Earth.

In keeping with this, I'll try to use coefficients alm that will give me a random field that is is isotropic on average. This seems hard at first, so let's just make a hand-waving argument. Looking at some pretty pictures of spherical harmonics, we can see that each spherical harmonic of degree l has about l stripes on it, independent of m. So let's try using alm's that depend only on l, and are constant if just m changes. [2]. Just for convenience, we'll pick this constant to be l to some power -p:
a_{lm} = l^{-p} \quad \textrm{ or}
h(\theta, \phi) = \sum_{l,m} N_{lm}(0,1) l^{-p} Y_l^m(\theta, \phi) \quad (2)

At this point I got bored & decided to see what a planet would look like if we didn't know what value of p to use. So below is a movie of a randomly generated "planet" with a fixed choice of random numbers, but with the power p changing.

As the movie starts (p=0), we see random uncorrelated heights on the surface [5]. As the movie continues and p increases, we see the surface smooth out rapidly. Eventually, after p=2 or so, the planet becomes very smooth and doesn't look at all like a planet. So the "correct" value for p is somewhere above 0 (too bumpy) and below 2 (too smooth). Can we use more observations about Earth to predict what a good value of p should be?

Observation 5: The elevation of the Earth's surface exists everywhere on Earth (duh).

So we're going to need our sum to exist. How the hell are we going to sum that series though! Not only is it random, but it also depends on where we are on the planet!

Rather than try to evaluate that sum everywhere on the sphere, I decided that it would be easiest to evaluate the sum at the "North Pole" at theta=0. Then, if we picked our coefficients right, this should be statistically the same as any other point on the planet.

Why do we want to look at theta = 0? Well, if we look back at the wikipedia entry for spherical harmonics, we see that
Y_l^m = \sqrt{ \frac{2 l + 1}{4\pi}\frac{(l-m)!}{(l+m)!}} e^{im\phi}P^m _ l(\cos\theta) \quad (3)

That doesn't look too helpful -- we've just picked up another special function Plm that we need to worry about. But there is a trick with these special functions Plm: at theta = 0, Plm is 0 if m isn't 0, and Pl0 is 1. So at theta = 0 this is simply:
Y_l^m(\theta = 0) = \bigg \{ ^{\sqrt{(2l+1)/4\pi},\, m=0}_{0,\,m \ne 0}

Now we just have, from every equation we've written down:
h(\theta = 0) = \sum_l \times l^{-p} \times \sqrt{(2l+1)/4\pi }\times N(0,1) \]\[
\quad \qquad =  \times \frac 1 {\sqrt{4\pi}} \times \sum_l N(0,l^{-2p}(2l+1)) \]\[
\quad \qquad =  \times \frac 1 {\sqrt{4\pi}} \times N(0,\sum_l l^{-2p}(2l+1) ) \] \[
\quad \qquad = \times \frac 1 {\sqrt{4\pi}} \sqrt{\sum_l l^{-2p}(2l+1)}
\times N(0,1) \] \[
\quad \qquad \sim \sqrt{\sum_l l^{-2p+1}} N(0,1) \qquad (4)

So for the surface of our imaginary planet to exist, we had better have that sum converge, or -2p+1 < -1 (p > 1).

And we've also learned something else!!! Our model always gives back a Gaussian height distribution on the surface. Changing the coefficients changes the variance of distribution of heights, but that's all it does to the distribution. Evidently if we want to get a non-Gaussian distribution of heights, we'd need to stretch our surface after evaluating the sum.

Well, what does the height distribution look like from my simulated planets? Just for the hell of it, I went ahead and generated ~400 independent surfaces at ~40 different values for the exponent p, looking at the first 22,499 terms in the series. From these surfaces I reconstructed the measured distributions; I've combined them into a movie which you can see below.


As you can see from the movie, the distributions look like Gaussians. The fits from Eq. (4) are overlaid in black dotted lines. (Since I can't sum an infinite number of spherical harmonics with a computer, I've plotted the fit I'd expect from just the terms I've summed. ) As you can see, they are all close to Gaussians. Not bad. Let's see what else we can get.

Observation 6: According to some famous people, the Earth's surface is probably a fractal whose coastlines are non-differentiable.

This means that we want a value of p that will make our surface rough enough so that its gradient doesn't exist (the derivative of the sum of Eq. (2) doesn't converge). At this point I'm getting bored with writing out calculations, so I'm just going to make some scaling arguments.

From Eq. (3), we know that each of the spherical harmonics Ylm is related to a polynomial of degree l in cos(theta). So if we take a derivative, I'd expect us to pick up another factor of l each time. Following through all the steps of Eq. (4) we find
\vec{\nabla}h \sim \sqrt{\sum_l l^{-2p+3}}\vec{N}(0,1) \quad ,
which converges for p > 2. So for our planet to be "fractal," we want 1<p<2. Looking at the first movie, this seems reasonable.

Observation 7: 70% of the Earth's surface is under water.

On Earth, we can think of the points underwater as all those points below a certain threshold height. So let's threshold the heights on our sphere. If we want 70% of our generated planet's surface to be under water, Eq (4) and the cumulative distribution function of a Gaussian distribution tells us that we want to pick a critical height H such that
\frac 1 2 \left[ 1 + \textrm{erf}(H/\sqrt{2\sigma^2}) \right] = 0.7 \quad \textrm{or}
\] \[
H = \sqrt{2\sigma^2}\textrm{erf}^{-1}(0.4)
\] \[
\sigma^2 = \frac 1 {4\pi} \sum_l l^{-2p}(2l+1) \quad (5)\, ,
where erf() is a special function called the error function, and erf-1 is its inverse. We can evaluate these numerically (or by using some dirty tricks if we're feeling especially masochistic).

So for our generated planet, let's call all the points with a height larger than H "land," and all the points with a height less than H "ocean." Here is what it looks like for a planet with p=0, p=1, and p=2, plotted with the same Sanson projection as before.
Top to bottom: p=0, p=1, and p=2. I've colored all the "water" (positions with heights < H as given in Eq. (5) ) blue and all the land (heights > H) green.

You can see that the the total amount of land area is roughly constant among the three images, but we haven't fixed how it's distributed. Looking at the map above for p=0, there are lots of small "islands" but no large contiguous land masses. For p=2, we see only one contiguous land mass (plus one 5-pixel island), and p=1 sits somewhere in between the two extremes. None of these look like the Earth, where there are a few large landmasses but many small islands. From our previous arguments, we'd expect something between p=1 and p=2 to look like the Earth, which is in line with the above picture. But how do we decide which value of p to use?

Observation 8: The Earth has 7 continents

This one is more vague than the others, but I think it's the coolest of all the arguments.

How do we compare our generated planets to the Earth? The Earth has 7 continents that comprise 4 different contiguous landmasses. In order, these are 1) Europe-Asia-Africa, 2) North- and South- America, 3) Antartica, and 4) Australia, with a 5th Greenland barely missing out. In terms of fractions of the Earth's surface, Google tells us that Australia covers 0.15% of the Earth's total surface area, and Greenland covers 0.04%. So let's define a "continent" as any contiguous landmass that accounts for more than 0.1% of the planet's total area. Then we can ask: What value of p gives us a planet with 4 continents?

I have no idea how to calculate exactly what that number would be from our model, but I can certainly measure it from the simulated results. I went ahead and counted the number of continents in the generated planets.

The results are plotted above. The solid red line is the median values of the number of continents, as measured over 400 distinct worlds at 40 different values of p. The red shaded region around it is the band containing the upper and lower quartiles of the number of continents. For comparison, in black on the right y-axis I've also plotted the log of the total number of landmasses at the resolution I've used.

The number of continents has a resonant value of p -- if p is too small, then there are many landmasses, but none are big enough to be continents. Conversely, if p is too large, then there is only one huge landmass. Somewhere in the middle, around p=0.5, there are about 20 continents, at least when only the first ~23000 terms in the series are summed.

Looking at the curve, we see that there are roughly two places where there are 4 continents in the world -- at p=0.1 and at p=1.3. Since p=0.1 doesn't converge, and since p=0.1 will have way too many landmasses, it looks like a generated Earth will look the best if we use a value of p=1.3

And that's it. For your viewing pleasure, here is a video of three of these planets below, complete with water, continents, and mountains. [4]


[1] Since I wanted a random surface, I wanted to make the mean of each coefficient 0. Otherwise we'd get a deterministic part of our surface heights. I picked a distribution that's symmetric about 0 because on Earth the bottom of the oceans seem roughly similar in terms of changes in elevation. I wanted to pick a stable distribution & independent coefficients because it makes the sums that come up easier to evalutate. Finally, I picked a Gaussian, as opposed to another stable distribution like a Lorentzian, because the tallest points on Earth are finite, and I wanted the variance of the planet's height to be defined.

[2] We could make this rigorous by showing that a rotated spherical harmonic is orthogonal to other spherical harmonics of a different degree l, but you don't want to see me try.

[3]Actually p=0 should correspond to completely uncorrelated delta-function noise. (You can convince yourself by looking at the spherical harmonic expansion for a delta-function.) The reason that the bumps have a finite width is that I only summed the first 22,499 terms in the series (l=150 and below). So the size of the bumps gives a rough idea of my resolution.

[4] For those of you keeping score at home, it took me more than 6 days to figure out how to make these planets.


  1. Random thoughts: few non-coastal plains and few coastal mountain ranges. Nice work!

    1. Glad somebody noticed that! I think that there being few non-coastal plains & few coastal mountains is a result of assuming uncorrelated coefficients. For all the ``planets'' I've looked at with different values of p, they all display the same problem/feature.

    2. I came across this article when I was searching for Spherical distribution for 3D montecarlo. But got entirely different but very interesting stuff in this. You did Awesome job, Bravo!

  2. good. thanks for the tutorial. it is very detailed and helpful

  3. Golf cart lift kits are methods to add many height along with improved suspension on your existing tennis cart. There are several different companies and sorts of golf carry lift products. In order to ascertain which golfing cart raise kit is befitting your type of golf basket, you ought to locate some form of chart as well as tool in which lists models of golf carts and the ggolf cart lift kits which have been compatible with each type.

  4. Great post, you have pointed out some excellent points, I as well believe this is a very superb website.
    Plastic Flow Meter

  5. ...interstellar travel constant acceleration (Sun-Deneb: 1000g)... Earth...the 2 ships that will go formation flying for mutual assisting if there are problems...indestructible structures made of Hexapentas material, awaiting in airport the arrival of passengers... Day 1: zero-speed... THE SHIPS TAKEOFF►... navigation computer places on screen the spacecraft in the center of sphere...spherical\tridimensional\spatial Heading: Deneb... Antimatter rocket engines...ON... Here we go...goooooo!...1g...10g...100g...constant acceleration cruise: 1000g (9.8 kms/sec²)... Inside of living areas (the same as going submerged in water: constant acceleration downwards...less...constant thrust, constant acceleration, from water upwards)...the gravitational transformers, perfectly synchronized with the acceleration, running: 1000g constant acceleration toward the floor ↓↓(motors)↓↓...less...999g constant acceleration toward the ceiling ↑↑(gravitational transformers)↑ = 1g constant acceleration toward the floor↓... 8.5 hours: light-speed = 1c...the fusion reactor as an artificial sun illuminating the immense Vital Support Gardens to lowering, from their comfortable apartments, cheerful passage to the pool...the electromagnetic shield anti-radiation...antigravity fields generator run forward, working: light objects away from the path of the ship...and trajectory ship away from the heavy objects...superluminal-speed > 1c... 42.5 hours: reaches hyperluminal-speed = 5c... Day 508: Half Journey...1000 light years...high hyperluminal-speed = 1435.39c... OFF engines...a few minutes of weightlessness during maneuver...the ship rotates 180º around its axis...motors ON again and... ◄STARTS TO BRAKE... Day 1017 (2.79 years): End Path party...2000 light The forever young passage of the 1st Immortal Generation (3D Bioprinting...Telomerase...modified Biological Timers...) disembarks at destination: an extra-stellar planet wich came errant to orbit of Deneb giant.

  6. ...interstellar travel (thousands G of constant acceleration)... gravitation, dimensions and inertia: the Matter... The electromagnetic radiation: quanta (pieces) of energy in discrete amounts (minimum and independents), photons, which move in waves... Why energy quanta (photons) are attracted by gravitational fields of "stellar lenses" curving their trajectories and move on)...go ahead because they have inertia...and are attracted because the Energy also has Gravitation, but does not manifest because it no emits Graviton, or emits in imperceptible degree, because it has to be condensed enough like in matter for that, but it is sensitive to Gravitation from Matter emits... Energy and Gravitation "associated"...and...the Gravitational Force is another manifestation, unknown yet, from Energy...(no "curvature" of that called "space-time" relativistic)... When was finished contraction, F=G((m1*m2)/d²), the pre-early Universe; (no relativistic "singularity", but it was a Sphere with dimensions and outside the infinite Empty Space); was need have as its starting point towards the expansion, one first action: when all Universe´s Matter was almost infinitely compressed and hot with all subatomic particles in direct contact without empty space among them, then →Gravity← becomes ←Antigravity→, those hyper-dense Matter becomes Energy ("infinites": Temperature + Energy + suddenly Antigravity = Big-Bang)...already again how a New Universe in expansion cooling, appears again the Matter from Energy...and again antigravity becomes normal Gravity... Its particle exchange, the theoretical Graviton yet...maybe "quanta of gravitational energy"... The Energy (e=mc²) is sensitive to Gravitation but hardly emits Graviton... The Matter. which is condensed Energy (m=e/c²), emits now Gravitons proportionally to its density... The Matter with normal density, naturally, emits few Gravitons (weak intensity of gravitational force). But...if someday could transform some of matter, not in Energy by annihilation with antimatter but entirely in Gravitational Force with an exponential Gravitons (or antigravitons) emission would have a way to amplify the Gravitational Force a certain amount of mass... Light (photons) can be directed wit a mirror...Gravity (gravitophotons) perhaps also can be directed someday the same such as light in a car´s headlight placed in ceiling looking downwards... This would be necessary for having artificial gravitational attraction only toward the ceiling in all the floors of the ship, each down´s floor must have less own gravity-power because receives gravitational attraction from all up´s floors. This can be difficult to observe in Universe because, although with unequal intensity, there are Gravitational Forces actuating from all directions. But in...gravitational transformers...for spacecrafts to thousands G of constant acceleration...

    1. (2)...interstellar travel (thousands G of constant acceleration)... The Electromagnetic and Gravitational Forces have some obvious similarities (fulfill the Law of Inverse of the Square) as both are "inversely proportional to the square of the distance that separate"...distance between the poles of two magnets the Electromagnetic...and distance between the centers of two masses the Gravitational... Just as now can manipulate the Electromagnetic Force transforming it into different voltages, currents and electrical powers...electrical transformers...perhaps someday when that mysterious Gravitational Force reveal its secrets, can also transform and amplify it (converting the mass exponentially in Gravitons) somehow in different intensities and gravitational powers...gravitational transformers... Spacecrafts with indestructible structures and the distribution of a building of housing flats (where are the foundations would the nozzles)... Antigravitational fields directed forward would go shoving and accelerating the space´s dust until reaching (the dust) almost ship-speed the dust´s cloud matter-feeder ("fuel") always placed at ship´s front, leading dust with Gravitational peripheral fields; already the dust with slow differential-speed between dust and ship (each dust´s particle already does not is a bullet); towards the ship converter matter-antimatter, where at half of matter would transform in antimatter to feed the insatiable demand for the voracious annihilation rocket engines ejecting the appreciable mass of the Pions, produced in large quantities, through the nozzles to the speed of light...a rocket´s thrust=linear momentum (mass*speed) ejection by nozzle in each second...and maintaining indefinitely to the ship at thousands G of constant acceleration... Powerful gravitational fields inside the living areas would counteract the HUGE CONSTANT ACCELERATION of the ship, keeping them to 1G...

  7. ...interstellar travel constant acceleration (times to reach light-speed)... 1g (9.8 mts/sec²)=354 days ...↓... 0.5g (4.9 mts/sec²)=1.94 years ...↓... 0.25g (2.45 mts/sec²)=3.88 years ...↓... 0.1g (98 cms/sec²)=9.7 years ...↓... 0.01g (9.8 cms/sec²)=97 years ...↓... 0.001g (9.8 mms/sec²)=970 years ______ Huge constant acceleration the ship, living areas to 1g: inside the living areas..."the same as going submerged in water"... ↓↓↓↓g (motors)...less...↑↑↑g (gravitational transformers)=↓g... ...↨... 2g (19.6 mts/sec²)=177 days ...↨... 10g (98 mts/sec²)=35.4 days ...↨... 100g (980 mts/sec²)=3.54 days ...↨... 1000g (9.8 kms/sec²)=8.5 hours ...↨... 10000g (98 kms/sec²)=51 minutes ...↨... 50000g (490 kms/sec²)=10.2 minutes ...↨... 100000g (980 kms/sec²)=5.1 minutes ...↨... 1 million G (9800 kms/sec²)=30.6 seconds... ((typewrite: interstellar travel constant acceleration))

  8. ...interstellar travel (laser)... L=2*3.14*r... that beam of a laser pointer that turn 180º, from horizon to horizon, in 1 second and 4 imaginary semicircular screens located, for example, to the distance from Earth... 1: 95493 kms... 2: Moon (384403 kms)... 3: Sun (150 Million kms, 1.5*10^8)... 4: Andromeda (2 million light years...19 Trillion kms, 1.9*10^19, Universal Large Scale: 1 billion=1 million of millions, 1*10¹²)... The mark of laser pointer to reach, would move alongside each screen to an angular speed of 180º/second and a linear speed, based on the light, of...(3.14*r)/c... km 95493= 1c... Moon= 4c... Sun= 1570c... Andromeda galaxy= 198 Billion*c (1.98*10^14)... If turn laser pointer to the right, why is going to bend each laser beam to the left, not to exceed c, from km 95493, radius "frontier" in which the speed of laser mark on the screen = c?... It could also be that a higher speed than speed of light, electromagnetic radiation becomes invisible and undetectable with current technology, entering a new and still unknown "superluminal dimension" wherein speed of light would be the minimum possible speed... Decreasing magnitudes have a close limit, but growing...have an "infinite" limit, to our knowledge, understood here as "infinite", at least, an exorbitant with the TEMPERATURE (-273 ºC→"infinite")... MATTER (quark→"infinite")... TIME (0→"infinite")... SPACE (0→"infinite")... why it would be SPEED different...(0 kms/sec→"infinite")... A rapidly rotating star, as a pulsar, is something similar as laser pointer example and the fact that it has not observed any radiation rate >c...seems to be relativistic reason. There is something that Science does not know...yet.