Least perimeter partition of the disc into N bubbles of two different areas

We present conjectured candidates for the least perimeter partition of a disc into \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$N \le 10$\end{document} connected regions which take one of two possible areas. We assume that the optimal partition is connected and enumerate all three-connected simple cubic graphs for each N. Candidate structures are obtained by assigning different areas to the regions: for even N there are N/2 bubbles of one area and N/2 bubbles of the other, and for odd N we consider both cases, i.e. in which the extra bubble takes either the larger or the smaller area. The perimeter of each candidate structure is found numerically for a few representative area ratios, and then the data is interpolated to give the conjectured least perimeter candidate for all possible area ratios. For each N we determine the ranges of area ratio for which each least perimeter candidate is optimal; at larger N these ranges are smaller, and there are more transitions from one optimal structure to another as the area ratio is varied. When the area ratio is significantly far from one, the least perimeter partitions tend to have a “mixed” configuration, in which bubbles of the same area are not adjacent to each other.


Introduction
Due to their structural stability and low material cost, energy-minimizing structures have a wide array of applications [1,2]. An example is the Beijing Aquatics Centre, which uses slices of the Weaire-Phelan foam structure [3] to create a lightweight and strong but beautiful piece of architecture.
The Weaire-Phelan foam is a solution to the celebrated Kelvin problem, which seeks the minimum surface area partition of space into bubbles of equal volume [4]. This builds upon the well-known isoperimetric problem concerning the least perimeter shape enclosing a given area [5]. Extending this idea to many bubbles with equal areas has led to further rigorous results for optimal structures, for example the proof of the honeycomb conjecture [6], the optimality of the standard triple bubble in the plane [7] and of the tetrahedral partition of the surface of the sphere into four regions [8].
If the areas of the bubbles are allowed to be unequal, then the problem of seeking the configuration of least perimeter is more difficult. For N = 2 regions in R 3 , the double bubble conjecture has been proved [9], and, in the plane, the extension of the honeycomb to two different areas (bidisperse) has led to conjectured solutions [10]. There has also been some experimental work that sought to correlate the frequency with which different configuraa e-mail: foams@aber.ac.uk tions of bidisperse bubble clusters (which, to a good approximation, minimize their surface area [1]) were found with the least perimeter configuration [11].
Minimal perimeter partitions of domains with a fixed boundary have also generated interest, for example a proof of the optimal partition of the disc into N = 3 regions of given areas [12], and many numerical conjectures, e.g. [13][14][15][16]. Such results may lead to further aesthetically pleasing structures like the Water Cube but that are truly foam-like, including their boundary, rather than being unphysical sections through a physical object.
It would be useful to seek patterns within the class of conjectured minimizers, so that we might be able to make general statements about the arrangement of regions within the domain. For example, do smaller bubbles tend to cluster together, since a small bubble adjacent to a large bubble might deform the larger bubble, increasing its perimeter? How does any such clustering depend on the differences in areas between bubbles? Fortes et al. [10,17] considered partitions of the plane into bubbles of two different areas and estimated the perimeter "cost" of placing bubbles of one area next to another. They predict that mixed small and large bubbles are favourable when the difference in areas is large [17] while when the two areas are similar the bubbles of each size tend to segregate, a so-called "sorted" configuration. Mechanical agitation often "shuffles" bidisperse foams towards such perimeterminimizing configurations [18][19][20].  In this work we seek to generate and test, in a systematic way, candidate partitions of domains with fixed boundary, and to examine the optimal patterns found to support conjectures about mixed and sorted configurations.
Due to the complexity, and in particular the large number of candidates, we restrict ourselves to a twodimensional (2D) problem. Thus, we enumerate all partitions of a disc and evaluate the perimeter of each one to determine the optimal configuration of the bubbles.
As the number of bubbles N increases then so does the complexity of the system and for N ≥ 5 numerical methods must be employed. For example, fig. 1 shows the two three-connected "simple" partitions of the disc into N = 5 bubbles with equal area. The difference in perimeter comes from the different structural arrangements of the arcs separating the bubbles. If we allow three bubbles to have one area and the other two a different area then there are 20 possible structures. When N = 10 this number increases to 314748.
We will use combinatorial arguments to enumerate the graphs corresponding to all possible structures. We recognise that all structures must obey Plateau's laws [21], a consequence of perimeter minimization [22], which state that edges have constant curvature and meet in threes at an angle of 2π/3. Rather than applying these directly, we will rely on standard numerical minimization software to determine the equilibrated configuration for each choice of N and areas.

Enumeration and evaluation of candidate structures
We consider each structure to be a simple, three-regular (cubic), three-connected planar graph ( fig. 2). This allows is to use graph-generating software to enumerate all the relevant partitions of the disc, and then to minimize each one's perimeter to determine the optimum. We therefore hypothesise a one-to-one correspondence between these graphs and the candidate solutions to the least perimeter partition. The assumption of planarity is natural, since these graphs must be embeddable in the 2D disc. The assumption that the graphs are three-regular follows from Plateau's laws. We assume that the graphs are simple and three-connected because any two edges sharing two vertices can be decomposed into a configuration with lower perimeter P . An example is shown in fig. 3: moving the lens-shaped bubble to the edge of the disc results in a change in topology and a reduction in perimeter. A similar reduction in perimeter can be achieved in structures with more bubbles by moving a lens towards a threefold vertex and performing the same change in topology.
Our conjectured minimizers also rely on the assumption that the bubbles in each candidate are connected. That is, no bubble may be split into more than one part in an attempt to reduce the total perimeter. This is something that is difficult to prove to be the case in general, but which has been shown to be true of all known minimizers [23].
We use the graph-enumeration software CaGe [24] to generate every graph (using the "3-regular Plane Graphs" generator with its default values), and an associated embedding, for each value of N (which corresponds to graphs with 2N − 2 vertices). This information is stored as a list of vertices, each with an (x, y) position and a list of neighbours. The number of graphs for each N is given in table 1.
The Surface Evolver [25] is a finite element software for the minimization of energy subject to constraints. We convert the CaGe output into a 2D Surface Evolver input file [15], in which each edge is represented as an arc of a circle (circular arc mode) and the relevant energy is the sum of edge lengths. The cluster is confined within a circular constraint with unit area, and we set a target area for each bubble. The Evolver's minimization routines are then used to find a minimum of the perimeter for each topology and target areas.
If an edge shrinks to zero length during the minimization, this is not a topology that will give rise to a stable candidate, since four-fold vertices are not minimizing. We therefore allow topological changes when an edge shrinks below a critical value l c (we use l c = 0.01, which is less than 1/50th of the disc radius). This prevents time-consuming calculation of non-optimal candidates, but does result in some solutions being found repeatedly as the result of different topological changes on different candidates. Fig. 3. The central two-sided bubble in this non-simple, two-connected structure (left) can be moved so that one of its vertices touches the boundary (middle) without changing the perimeter P = 5.665 of this configuration. Once there, a change in topology results in a drop in the perimeter to P = 5.196 and a simple, three-connected, state (right). Table 1. For each number of bubbles N we show the number of simple, cubic, three-connected graphs, the number of permutations of the two possible areas (for odd N this is half of the number of structures tested), and then the product, which is the number of candidates whose perimeter we evaluate. The last five columns give the number of distinct realizable structures found after minimization, for each area ratio.  Our aim is to consider bidisperse structures, in which each bubble can take one of two possible areas. We define the area ratio A r to be the ratio of the area of the large bubbles to the area of the small bubbles, so that A r > 1. When the smaller bubbles are very small, the precise area ratio changes the total energy only very little, so we consider A r up to 10. (The highest area ratio at which we find a change in the topology of the optimal structure is A r = 8.35.) To reduce the number of possible candidates, we stipulate that the number of bubbles of each area are equal (when N is even) or (when N is odd) as close as possible. For example, for N = 3, there is only one possible graph ( fig. 2), in which three lines meet together in an internal vertex, as for the monodisperse case. Since N is odd we consider 3 21 and 3 12 separately. In the first case there are three possible permutations of the areas assigned to the three bubbles, but all three are clearly equivalent through a rotation, so there is only one candidate for which the perimeter must be evaluated. In the second case there are also three possible permutations of the area, but again only one candidate needs to be minimized.
The number of graphs and the number of area permutations rises rapidly. We therefore treat only values of N between 4 and 10. The number of candidates that we evaluate and the number of structures that are actually realized is shown in table 1.

Least perimeter candidates at representative area ratios
The perimeter P decreases quite strongly with increasing area ratio, because small enough bubbles make only a small perturbation to a structure with lower N , and structures with lower N have lower P . Although the average area of each bubble is fixed (at 1/N ), the polydispersity increases with A r . A general measure of polydispersity for bubbles i with areas A i N is where denotes an average over all bubbles i. Note that with this definition p = 0 for a monodisperse partition. For a partition with N L large bubbles this becomes We expect the perimeter to decrease as 1/(1 + p) [26], and so to help distinguish different candidates for given N over a range of area ratio A r , we plot P (1 + p) in the following. Figures 5-11, for N = 4 to 10 respectively, show the scaled perimeter of the structures analysed. The optimal perimeter for each N and each A r is highlighted with a thick line, the transitions between structures are indicated, and the least perimeter structures themselves are shown according to the area ratio at which they are found.
We start by investigating area ratios A r = 2, 4, 6, 8 and 10. For N = 4 and 5 there is no change in the topology of our conjectured least perimeter structure as the area ratio changes; see figs. 5 and 6. For N = 4 the two smaller bubbles never touch, and lie at opposite ends of a straight central edge. For N = 5, for both possible distributions of large and small bubbles, the optimal pattern always consists of two three-sided bubbles whose internal vertices are connected to the other internal vertex, which itself has one other connection to the boundary of the disc. That is, in neither case does the optimal candidate have an internal bubble.
For N ≥ 6 there are transitions between different structures as the area ratio A r changes. We therefore interpolate between these values of area ratio to determine the critical values of A r at which the changes in topology of the least perimeter candidate occur for each N .
We do this by taking each of the structures that was found for each area ratio A r = 2, 4, . . . and change the area ratio in small steps (of 0.05). For each of these candidates we find and record the perimeter. (For N = 10 we do this only for the fifty or so best candidates for each value of A r , since there are so many candidates which are far from optimal for any area ratio.) For candidates whose initial area ratio was 2, 4 or 6 we decreased the area ratio to 1.1 and then increased it up to 10. For candidates whose initial area ratio was 8 or 10 we increased the area ratio up to 10 before slowly decreasing it down to 1.1. We are therefore able to confirm that at low enough area ratio, i.e. close to 1, we recover the optimal structures found in the monodisperse case [27].
This procedure generates a few extra optimal structures that are missed by the first sampling of the area ratios, for example between A r = 2 and 4 for N = 8, 9 54 For N = 6 we find that the topology of the least perimeter candidates with area ratios A r ≥ 2.6 are the same (fig. 7). For area ratios less than this value the topology is that of the optimal candidate in the monodisperse case [27].
The two different cases for N = 7 behave differently ( fig. 8) Fig. 11. The perimeter P of the least perimeter candidates for N = 10 at different area ratios. We only show the perimeter corresponding to the fifty best candidates for each area ratio. ent topologies, with a transition to a new candidate at a surprisingly high area ratio of 8.4.
The least perimeter structure with N = 8 bubbles has the monodisperse topology for A r < 2.6 ( fig. 9); there is one further transition at A r = 3.9, giving three different optimal structures. For N = 9 the results are richer ( fig. 10), in the sense that the system explores more possible states as the area ratio changes. For 9 54 we find five different topologies, while for 9 45 there are four. In the latter case the structure found for A r = 2 is different to the monodisperse one [27], and there is a transition to that structure at a low area ratio around 1.8.
Finally, for N = 10 ( fig. 11) we again find a candidate for A r = 2 that differs from the monodisperse case and a transition at even lower area ratio. In total there are five different topologies.

Analysis of patterns
The critical area ratios at which there is a transition between optimal structures are summarised in fig. 12. Most are found at intermediate values of the area ratio, roughly between A r = 2.5 and 4, although this broadens slightly with increasing N . There is also a single point at high A r , for N = 7, which corresponds to moving a small bubble from the boundary of the disc to the centre, and hence to a symmetric state. It is perhaps surprising that this highly symmetric state is not optimal at lower area ratio, since many of the least perimeter structures are symmetric.
The images in figs. 5-11 also hint at an evolution from the small bubbles clustering together at low area ratio to being separated from each other by the large bubbles at high area ratio. This is exemplified by the highly symmetric structures for 7 34 and 9 45 at high area ratio. We quantify this observation by counting the proportion of edges E LS separating large from small bubbles in each least perimeter structure. A structure with a higher value of E LS has less clustering [18]. The data in fig. 13 bears out this observation: for N ≥ 6 and small area ratio the  Fig. 13. The proportion of edges separating large from small bubbles E LS in the least perimeter candidate for each N . The annulus size is proportional to the area ratio, and the thickness of each annulus represents the range of area ratios for which the candidate with the corresponding value of E LS is conjectured to be optimal. (For example, for N = 6 we plot an annulus with (in arbitrary units) inner radius 1 and outer radius 2.6 (from fig. 12) at E LS = 0.5, and another with inner radius 2.6 and outer radius 10 at E LS = 0.66.) Data for the cases of an extra large or an extra small bubble (N = 7, 9) is displaced slightly to reduce overlap. (Note that two adjacent structures, for N = 9 and 10, have the same value of E LS , resulting in two pairs of the smallest annuli being combined.) value of E LS is lower than for large area ratio. Moreover, a horizontal line at E LS ≈ 0.6 separates all configurations with area ratio less than about 4 from those (above the line) with larger area ratio. Within each of these two parts of the graph, it is not possible to distinguish a trend.

Conclusions
We have enumerated all candidate partitions of the disc with N ≤ 10 bubbles with one of two different areas, and determined, for each area ratio, the partition with least perimeter. The results show an increasing number of transitions between the different optimal structures found for varying area ratio as N increases, mostly at low area ratio. Further, in the least perimeter partitions at small area ratio the smaller bubbles are clustered together, while at large area ratio the small bubbles are separated by large bubbles. Transitions between such mixed and sorted configurations often occur as a consequence of some agitation [20].
The procedure described here should translate directly to least perimeter partitions of the surface of a sphere, since to enumerate candidates to that problem we are able to use the same graphs and consider the periphery of the graph to form the boundary of one further region. Thus the candidates for the disc with N bubbles are also the candidates for the sphere with N + 1 regions. In general, our preliminary results indicate, as for the monodisperse case [15], that the least perimeter arrangement of regions on the sphere is different to the corresponding optimal partition of the disc.