Updated 15 January 2023
How would bubbles pack together, to give the least possible amount of surface film between them? This question has not been answered yet.
The problem is actually a little more precise than that, as the obvious solution is one giant spherical bubble, but that is clearly not what we have in mind. What about bubbles of a given size? Well, why not allow bubbles of different sizes if that turns out to be more efficient? What we really mean is bubbles at a given density, where the bubble density is the number of bubbles per unit volume.
The faces and edges of all bubbles have to meet at the correct angles to minimise the surface energy, as described by Plateau's laws. For the same reason, bubbles always meet three to an edge and four to a vertex. The vertex angle is the same for all bubbles, and is known as the tetrahdral angle. However polyhedra with only this face angle do not exist, so some of the faces and edges must always curve a little to reach the correct angle. The faces are minimal surfaces stretching between the edges.
In the nineteenth century Lord Kelvin considered the simplest case of identically sized and shaped bubbles. He discovered a packing of 14sided bubbles or tetrakaidecahedra, based on the truncated octahedron but with slightly curved faces. Kelvin conjectured that this was the most efficient solution, but he could neither prove nor disprove his conjecture.
Vertex coordinates and a printable net of the truncated octahedron are given on the right. Both give the same bubble volume as the WeairePhelan bubbles further down.
Truncated octahedron x y z 3.14980 0 6.29961 0 3.14980 6.29961 3.14980 0 6.29961 0 3.14980 6.29961 6.29961 0 3.14980 0 6.29961 3.14980 6.29961 0 3.14980 0 6.29961 3.14980 6.29961 3.14980 0 3.14980 6.29961 0 3.14980 6.29961 0 6.29961 3.14980 0 6.29961 3.14980 0 3.14980 6.29961 0 3.14980 6.29961 0 6.29961 3.14980 0 6.29961 0 3.14980 0 6.29961 3.14980 6.29961 0 3.14980 0 6.29961 3.14980 3.14980 0 6.29961 0 3.14980 6.29961 3.14980 0 6.29961 0 3.14980 6.29961 

A hundred years went by before Denis Weaire and Robert Phelan disproved the Kelvin conjecture by finding a better solution. They asked the slightly less strict question, how would equal sized, but not necessarily identically shaped, bubbles pack? The best answer they could find was a mixture of 12 and 14sided bubbles, also having slightly curved surfaces.[3]
In the WeairePhelan packing, two irregular pentagonal dodecahedra (12sided) and six tetrakaidecahedra (of a rather different shape from Kelvin's, which also has 14 faces) form a translation unit with a lattice periodicity which is simple cubic. The illustration shows the dodecahedra as wire frames and the tetrakaidecahedra as solid. The dodecahedra do not touch each other, but are entirely surrounded by tetrakaidecahedra.
The dodecahedron has pentagonal faces with sides of unequal length, so it is "pentagonal" to distinguish it from the rhombic dodecahedron and "irregular" to distinguish it from the regular (Platonic) dodecahedron.
Another way to visualise the packing is to note that the hexagonal faces of the tetrakaidecahedra are truly flat and meet to stack them in long rods. A set of rods lies parallel to each of the three orthogonal axes of space. These three sets of rods interlace, touching on their eight smaller pentagonal faces but leaving dodecahedral voids between the larger pentagons.
More about bubble packings can be found in [4].
Since then no better solution has been found, though at least one has come very close. Real WeairePhelan bubbles also resisted experimental demonstration for almost twenty years, but a foam was eventually made.[1] Proof of the optimal solution remains as elusive as ever, we still do not know if this is it.
We can flatten out all the bubble faces to form conventional polyhedra. These faces lie orthogonal to the lines joining adjacent cell centroids, and to retain equal cell volumes some vertices must move very slightly. The geometry of the resulting spacefilling is in all other respects identical to Weaire and Phelan's: the differences are, in any case, barely noticeable.
The tetrakaidecahedron may now be seen an example of the Goldberg polyhedron, which helps to distinguish it from the one associated with Kelvin, the truncated octahedron. The dodecahedron has the same symmetry seen in crystals of iron pyrites and so is sometimes called a pyritohedron, although its angles are different.
Vertex coordinates and printable nets for card models of the two polyhedra are given below. The coordinate scale gives a lattice period, i.e. the size of the repeating pattern, of 20.
Dodecahedron x y z 3.14980 0 6.29961 3.14980 0 6.29961 4.19974 4.19974 4.19974 0 6.29961 3.14980 4.19974 4.19974 4.19974 4.19974 4.19974 4.19974 0 6.29961 3.14980 4.19974 4.19974 4.19974 6.29961 3.14980 0 6.29961 3.14980 0 6.29961 3.14980 0 6.29961 3.14980 0 4.19974 4.19974 4.19974 0 6.29961 3.14980 4.19974 4.19974 4.19974 4.19974 4.19974 4.19974 0 6.29961 3.14980 4.19974 4.19974 4.19974 3.14980 0 6.29961 3.14980 0 6.29961 
Tetrakaidecahedron x y z 3.14980 3.70039 5 3.14980 3.70039 5 5 0 5 3.14980 3.70039 5 3.14980 3.70039 5 5 0 5 4.19974 5.80026 0.80026 4.19974 5.80026 0.80026 6.85020 0 1.29961 4.19974 5.80026 0.80026 4.19974 5.80026 0.80026 6.85020 0 1.29961 5.80026 4.19974 0.80026 0 6.85020 1.29961 5.80026 4.19974 0.80026 5.80026 4.19974 0.80026 0 6.85020 1.29961 5.80026 4.19974 0.80026 3.70039 3.14980 5 0 5 5 3.70039 3.14980 5 3.70039 3.14980 5 0 5 5 3.70039 3.14980 5 
The polyhedral spacefilling is related to what crystallographers call an A15 FrankKasper phase. These phases are characterised by the arrangement of atoms (here represented by the cell centroids) at the vertices of planar arrays of triangles and hexagons. A natural example may be found in certain gas hydrates, in which the gas molecules lie at the centres of the cells while the oxygen atoms of water form the vertices of the lattice and the hydrogen atoms lie along the edges. It is described by Williams as the 12 Å cubic gas hydrate, where the angstrom Å is a unit of length equal to 10^{−10} m.[4]
The WeairePhelan polyhedra approximate the Voronoi cells of such a crystal. If you place a small sphere at each trapped gas molecule in the array and expand them all together until they meet to form a boundary, these boundaries will be flat faces. If you keep expanding the spheres until no intervening space is left, the flat faces all join up and the spheres become polyhedra known as Voronoi cells.However the classic Voroni construction creates faces at equal distances frrom the adjacent cell centroids, leading in this case to 12hedra which have a slightly different volume from the 14hedra. This can be adjusted by weighting the growth rates of the cells slightly, analogous to adjusting the bubbles' internal pressures. The equalvolume polyhedral spacefilling is therefore a weighted Voroni construction.[2]
The differences between the bubbles and the two polyhedral constructions are small and barely noticeable. To be honest, updating this page many years after the event, I can no longer recall which polyhedral construction I use here.
Since it has flat faces, the polyheral spacefilling is better suited to making card models than the actual bubbles are. Downloadable nets are provided above.
The polyhedra are also easier to design for 3D printing than the bubbles! Tip: it can be worth making them hollow with a fairly thin wall, so as to use less material. This saves cost and time, especially if you are printing off large numbers, and also makes them lighter.