Mapping the Fibonacci lattice (aka Golden Spiral, aka Fibonacci Sphere) onto the surface of a sphere is an extremely fast and effective approximate method to evenly distribute points on a sphere. I show how small modifications to the canonical implementation can result in notable improvements for nearest-neighbor measures.

First published: 07 June, 2020

Last modified: 09 June, 2020

This is an updated version of my most popular post “Evenly distributing points on a sphere”.

# 0. The Fibonacci lattice

The problem of how to evenly distribute points on a sphere has a very long history. It is of huge signficance in many areas of mathematics, physics, chemistry and computing – including numerical analysis, approximation theory, crystallography, viral morphology, electrostatics, coding theory, and computer graphics, to name just a few.

Unfortunately, except for a small handful of cases, it still has not been exactly solved. Therefore, in nearly all situations, we can merely hope to find near-optimal solutions to this problem.

Of these near-optimal solutions, one of the most common simple methods is one based on the Fibonacci lattice or the Golden Spiral. This method was apparently inspired by phyllotaxis (which describes the seed distribution on the head of a sunflower or a pine cone).

Furthermore, unlike many of the other iterative or randomised methods such a simulated annealing, the Fibonacci spiral is one of the few direct construction methods that works for arbitrary $n$.

The usual modern definition of the Fibonacci lattice (see figure 2, top row), which evenly distributes $n$ points inside a unit square $[0, 1)^2$ is:

$$ t_i = (x_i,y_i) = \left( \frac{i }{\phi} \% 1 , \frac{i}{n}, \right) \quad \textrm{for }\; 0 \leq i < n \tag{1}$$

where

$$\phi = \frac{1+\sqrt{5}}{2} = \lim_{n \rightarrow \infty} \left( \frac{F_{n+1}}{F_n} \right)$$

and the $(\cdot) \%1$ operator denotes the fractional part of the argument.

(This method also works very well for rectangles, too.)

To produce the Fibonacci spiral, aka the Golden Spiral (see figure 2, bottom row), you simply map this point distribution to a unit disk, via the equal-area transformation

$$ (x,y) \rightarrow (\theta, r) : \quad (2\pi x, \sqrt{y}) $$.

Similarly, these point sets can be mapped from the unit square $[0, 1)^2$ to the sphere $S^2$ by the cylindrical equal-area projection (see figure 2):

$$ (x,y) \rightarrow (\theta, \phi) : \quad \left(2 \pi x, \arccos(1-2y), \right) $$

$$ (\theta,\phi) \rightarrow (x,y,z) : \quad \left (\cos\theta \sin\phi, \sin \theta \sin \phi, \cos \phi \right) $$

Note that there are numerous conventions for labelling angles in spherical geometry. The above formula uses the notation that $\theta \in [0,2\pi]$ is the longitude and $\phi \in [0,\pi]$ is the angle from North pole.

That is, the Fibonacci lattice is a simple way to very evenly distribute points in a rectangle, disk or the surface of a sphere.

(See footnote #2, on why the Fibonacci lattice does not wrap on all edges, and why this is a good thing.)

Here’s a basic implementation of this in Python.

```
from numpy import arange, pi, sin, cos, arccos
n = 50
goldenRatio = (1 + 5**0.5)/2
i = arange(0, n)
theta = 2 *pi * i / goldenRatio
phi = arccos(1 - 2*(i+0.5)/n)
x, y, z = cos(theta) * sin(phi), sin(theta) * sin(phi), cos(phi);
```

This spiral and spherical mapping both suffer from two distinct but inter-related problems.

The first is that this mapping is area-preserving, not distance-preserving. If you are attempting to optimize distance measures such as those related to nearest-neighbor, then it is not guaranteed that such distance-based constraints and relationships will hold after an area-based projection.

The second, and from a pragmatic point of view possibly the trickiest to resolve, is that the spiral mapping has a singularity point at the centre (and the spherical mapping has a singularity at each pole). In the spherical case, consider two points very close to the pole but 180 degrees different in longitude. On the unit square, (and also on the cylindrical projection) they would correspond to two points that are near the top edge, but still quite distant from each other. However, when mapped to the surface of the sphere they could be joined be a very small great arc going over the north pole. It is this particular issue that makes many of the spiral mappings sub-optimal. This issue does not exist for the fibonacci lattice in $[0,1)^2$.

# 1. Optimizing minimum nearest-neighbor distance

There are many equally valid criteria that are used to describe how evenly distributed the points on the surface of a sphere are. These include:

- Packing and covering
- Convex hulls, Voronoi cells and Delaunay triangles,
- Riesz $s$-energy kernels
- Cubature and Determinants

It is critical to realise that an optimal solution based on one criteria is often not an optimal point distribution for another.

In this post, we only explore two of these measures. This main section looks at the most famous measure: *minimum* nearest-neighbor distance; while the next section briefly describes *average* nearest-neighbor distance. (See foonote #3, on exploring how to optimize the Fibonacci lattice for other measures).

The problem of finding the maximum nearest-neighbor is often referred to as “*Tammes’s problem*” due to the botanist Tammes, who searched for an explanation of the surface structure of pollen grains. It’s also often called Fejes Tóth’s problem. And the configuration corresponding to this optimal nearest neighbor is often called a spherical code.

**Platonic Solids**

One myth that is important to dispel is the belief that all the Platonic solids are the optimal configurations for their respectively corresponding values of $n$.

Unfortunately, this is not true. This is succintly described on the Wolfram mathworld site:

“For two points, the points should be at opposite ends of a diameter. For four points, they should be placed at the polyhedron vertices of an inscribed regular tetrahedron. There is no unique best solution for five points since the distance cannot be reduced below that for six points. For six points, they should be placed at the polyhedron vertices of an inscribed regular octahedron. For seven points, the best solution is four equilateral spherical triangles with angles of 80 degress. For eight points, the best dispersal is not the polyhedron vertices of the inscribed cube, but of a square antiprism with equal polyhedron edges. The solution for nine points is eight equilateral spherical triangles with angles of $\arccos(\frac{1}{4})$. For 12 points, the solution is an inscribed regular icosahedron.”

Let’s repeat this point, as it is not obvious and not well-known:

Two of the platonic solids, the cube ($n=8$) and the dodecahedron $(n=20)$ are **not** the optimal point configurations for minimum nearest-neighbor distance. And actually, they are far from it.

So before we explore how to improve the canonical Fibonacci lattice for arbitrary $n$, let’s set up some notation.

If, for each of the $n$ points, we define $\delta_i$ to be distance from point $i$ to its closest neighbor $(0 \leq i < n$). Then Tammes’s problem asks us to maximize the minimum of $d_i$ over all $n$ points. That is,

$$ \delta_{\rm{min}}(n) = \min_{0\leq i < n} \delta_i $$

This value decreases at a rate ~$1/\sqrt{n}$, so it is also useful to define the normalized distance, and also the asymptotic limit of the normalized distance as

$$ \delta^*(n) = \sqrt{n} \delta(n) ,\quad \quad \delta^* = \lim_{n \rightarrow \infty} \delta^*(n) $$

**Canonical lattice**

The standard implementation (as listed in the code above), divides the unit interval $[0,1)$ into $n$ equal subintervals, and then places a point at the midpoint of each of these intervals (aka the midpoint rule in gausssian quadrature). That is,

$$ t_i (\theta, \phi) = \left( \frac{i}{\phi}, \frac{i{+0.5}}{n} \right) \quad \textrm{for }\; 0 \leq i < n \tag{1}$$

This method produces a consistent value of $\delta^*_n=3.09$ for all $n$.

**Offset Fibonacci Lattice.**

A key insight to improving on Equation 1, is to realize that the $\delta_\rm{min}(n)$ always corresponds to the distance between the points $t_0$ and $t_3$, which are located near the poles. Thus, to improve $\delta(n)$ the points near the poles should be positioned farther apart.

To achieve this, we need to move (offset) all the points slightly farther away from the poles. This of course means, that almost all of them become slightly closer together. However, as they are not the bottleneck for optimizing $\delta_\rm{min}(n)$, then in this situation, this is tolerable.

Specifically, we can achieve this by defining a sub-interval of length $\epsilon$ on each end of the interval $(0,1)$ where there are no points, and then we equally space the $n$ points on the remaining interval $(\epsilon, 1{-}\epsilon)$ ensuring that there is a point at each end of this interval.

That is,

$$ t_i(\varepsilon) = \left( \frac{i{+ \varepsilon}}{n{-1}{+2\varepsilon}}, \frac{i}{\phi} \right) \quad \textrm{for }\; 0 \leq i < n \tag{2}$$

Note that $\varepsilon=\frac{1}{2}$ corresponds to canonical lattice, and $\varepsilon>\frac{1}{2}$ represents an increased gap near the poles.

Through methods such as Monte Carlo, for each $n$ we can empirically find the optimal value of $\varepsilon$ that optimizes packing distance. Let’s denote this optimal value as $\varepsilon^*(n)$.

The bad news is that there does not seem to be a simple closed expression for $\varepsilon^*(n)$. However, the good news is that some very elegant properties can be found.

The first is that the function $\varepsilon^*(n)$ is not a continuous function, but more like a step function, (significantly) increasing only at certain values of $n$. Furthermore, the corresponding values of $\varepsilon(n)$ are nearly always very close to a round number. More specfically, empirical analysis shows that the following table is an excellent description for the best values of $\varepsilon$:

\begin{array}{|rcccl|r|} \hline

&&&&& \varepsilon(n) \\ \hline

\quad & &n &<&24 & 0.33\\

24 &\leq &n &<&177 & 1.33\\

177 &\leq &n &< &890 & 3.33 \\

890 &\leq &n &<&11,000 & 10 \\

11,000 &\leq &n &<&39,000 & 27 \\

400,000 &\leq &n &<&600,000 & 75\\

600,000 &\leq &n & \quad && 214\\ \hline

\end{array}

A basic implementation of this offset-lattice is shown.

```
from numpy import arange, pi, sin, cos, arccos
n = 50
if n >= 600000:
epsilon = 214
elif n>= 400000:
epsilon = 75
elif n>= 11000:
epsilon = 27
elif n>= 890:
epsilon = 10
elif n>= 177:
epsilon = 3.33
elif n>= 24:
epsilon = 1.33
else:
epsilon = 0.33
goldenRatio = (1 + 5**0.5)/2
i = arange(0, n)
theta = 2 *pi * i / goldenRatio
phi = arccos(1 - 2*(i+epsilon)/(n-1+2*epsilon))
x, y, z = cos(theta) * sin(phi), sin(theta) * sin(phi), cos(phi);
```

Offsetting the points of the Fibonacci lattice slightly away from the poles (as defined in equation 2) produces a packing $\delta^*_{\rm{min}}$, that is up to 8.3% tighter than the canonical Fibonacci lattice.

\begin{array}{|r|rr|r|} \hline

n & \rm{canonical} & \rm{new} & \rm{increase} \\ \hline

10 & 3.09 & 3.13 & 0.6\%\\

20 & 3.09 & 3.12 & 1.1\%\\

50 & 3.09 & 3.21 & 4.1\%\\ \hline

100 & 3.09 & 3.25 & 5.1\%\\

200 & 3.09 & 3.27 & 5.7\%\\

500 & 3.09 & 3.30 & 6.6\% \\ \hline

1,000 & 3.09 & 3.31 & 7.0\%\\

2,000 & 3.09 & 3.33 & 7.6\%\\

5,000 & 3.09 & 3.34 & 7.8\% \\ \hline

10,000 & 3.09 & 3.34 & 8.0\%\\

20,000 & 3.09 & 3.34 & 8.0\%\\

50,000 & 3.09 & 3.35 & 8.2\%\\ \hline

\end{array}

(See footnote #4 on additional ways to modify the canonical lattice that increases $\delta_{\rm{min}}$ even further.)

**Comparison**

We can compare this result with some of the state-of-the-art methods which are typically complex and require recursive and/or dynamic programming are as follows (the higher the better).

\begin{array}{|lr|} \hline

\text{} & d^* \\ \hline

\text{Coulomb} & 3.37 \\

\text{Log Energy} & 3.37\\

\text{Offset Fibonacci Lattice} & 3.35 \\

\text{Zonal Equal Area} & 3.32 \\

\text{Max Determinant} & 3.19 \\

\text{Canonical Fibonacci Lattice} & 3.09\\ \hline

\end{array}

# 2. Corollary. Preventing large gaps

For large value of $n$, (especially $n>1000$), when optimzing packing distance, the optimal value of $\varepsilon$ is quite large, and the effect of this is a noticeable gap or void at each pole.

A simple technique to avoid this gap, is to place a point at each pole, and then place the remaining $n{-}2$ points as per equation 2 with an revised offset $\varepsilon_2$as per the table below.

Not only does this minimize the gap at each pole, but it also offers a tiny improvement on the offset lattice as defined in equation 2. The size of the improvement is very small – only up to 1%, and typically less than 0.1%.

Therefore, in regards to packing distance it probably is probably more for theoretical interest than for any practical purpose. However, in most cases it is probably worth it to ensure wide applicability for the resulting configuration.

\begin{array}{|rcccl|r|} \hline

&&&&& \varepsilon_2(n) \\ \hline

\quad & &n &<&80 & 2.66\\

80 &\leq &n &<&1000 &3.33\\

1000&\leq &n &< &40000& 10\\

40,000 &\leq &n & \quad && 25\\ \hline

\end{array}

# 3. Optimizing average nearest-neighbor distance

Although maximizing minimum nearest neigbor is the most studied and is extremly useful from a theoretical perspective, it might not be the most useful one for many engineering and computational applications.

Unfortunately, optimzing for minimum nearest neighbor, produces less than optimal configurations for **average** nearest neighbor. This is not surprising as the offset pushes all but a small handful of points away from the poles, and thus closer together.

It turns out that for all $n$, there is essentially a single value of $\varepsilon$ that offers the best configuration. Namely, $\varepsilon = 0.36$.

Note that this value is less than 0.5, so the points are actually pulled *closer* to the poles compared to the canonical implementation!

Although offsetting the points with $\varepsilon= 0.36$ produces suboptimal configurations for packing distance, it produces a optimal configurations for **average** nearest neighbor measure. This makes it ideal for most engineering and computing applications.

# 5. Summary

Mapping the Fibonacci lattice (aka Golden Spiral, aka Fibonacci Sphere) onto the surface of a sphere is an extremely fast and effective approximate method to evenly distribute points on a sphere.

However, by offsetting the points of the Fibonacci lattice slightly away from the poles (as defined in equation 2), results in a point onfiguration that is up to 8.3% tighter than the canonical Fibonacci lattice (in terms of packing distance).

\begin{array}{|r|rr|r|} \hline

n & \rm{canonical} & \rm{new} & \rm{increase} \\ \hline

100 & 3.09 & 3.25 & 5.1\%\\

200 & 3.09 & 3.27 & 5.7\%\\

500 & 3.09 & 3.30 & 6.6\% \\ \hline

1,000 & 3.09 & 3.31 & 7.0\%\\

2,000 & 3.09 & 3.33 & 7.6\%\\

5,000 & 3.09 & 3.34 & 7.8\% \\ \hline

10,000 & 3.09 & 3.34 & 8.0\%\\

20,000 & 3.09 & 3.34 & 8.0\%\\

50,000 & 3.09 & 3.35 & 8.2\%\\ \hline

\end{array}

For $n>100$, an improvement can be made beyond this, by initially placing a point at each pole, and then placing the remaining $n-2$ points as per equation 2 with an offset as per table in section 4. This not only (very sightly) improves minimal nearest packing, but it also prevents a large gap at each pole.

**Footnotes**

#1. Technically, a sphere only consists of all points that are equidistant from the centre. This set of all such points forms a surface. Therefore the phrase ‘the surface of the sphere’ is a tautology. However, in this post, I frequently use this phrase to ensure there is no confusion amongst those readers, regardless of their technical background.

#2. The $\%1%$operator in the Fibonacci lattice (see figure 1, top row) indicates that we are really mapping the unit square to a cylinder. It is often considered that this new space is a torus, but careful observation shows that points near the top edge and the bottom edge do not precisely align. In our case where we are mapping the lattice to the a disc or to the surface of a sphere, this actually works in our favor as it reduces the constraints and therefore allows a higher degree of freedom when finding good point separations.

This is one reason why mapping Fibonacci lattice to the surface of the sphere nearly always produces better results than mapping a typical two-dimensional low discrepancy sequence (where all the edges wrap) to the surface of the sphere.

#3. This post explores measures relating to nearest neighbor. To optimize the point configurations such that the volume of the convex hull is maximized, see my previous post, “Evenly distributing points on a sphere”.)

#4. Moving the points farther and farther away from the points, suggests that there is a threshold where it is constructive to place a point at each of the poles. Details of this were explored in the previous post on this topic “evenly distributing points on a sphere“. However, with these newly found optimal values of $\varepsilon$, this additional complexity only improves $\delta_{\rm{min}}$ by up to 0.2%, so for most applications it is probably not worth implementing.

**Key References**

- Saff: A Comparison of poopular point confugrations on $S^2$“:
- Evenly_distributed_points_on_sphere
- Saff-Kuijlaars; MathIntel97
- Project Euclid; euclid.em/1067634731
- Vanderhaeghe; Spherical Fibonacci Mapping
- Leopardi; Sphere talk
- Leopardi; Sphere-PhD-Thesis
- mathworld.(wolfram); SphericalCode
- Neil Sloane: Spherical Packings

***

My name is Dr Martin Roberts, and I’m a freelance Principal Data Science consultant, who loves working at the intersection of maths and computing.

“I transform and modernize organizations through innovative data strategies solutions.”

You can contact me through any of these channels.LinkedIn: https://www.linkedin.com/in/martinroberts/

Twitter:

@TechSparxhttps://twitter.com/TechSparxemail: Martin (at) RobertsAnalytics (dot) com

More details about me can be found

here.

Other Posts you might like: