# Planar network statistics for two-dimensional rupturing foams

Joseph Klobusicky<EMAIL_ADDRESS>Elif Onat  Department of
Mathematics, The University of Scranton, Scranton, PA, 18510  Vasilios
Konstantinou  Lam Research, Tualatin, OR 97062

###### Abstract

We conduct experiments on a class of two-dimensional foams generated through
compressing a three-dimensional soap foam between two glass plates. To induce
a spatially uniform rupturing process on foam boundaries, an additional plate
is heated and placed on top of the unheated plates. For 30 separate foam
samples, we record network statistics related to cell side numbers and areas
as the foam coarsens over a half-minute. We find that the Aboav law and a
quadratic Lewis Law, two commonly used relations between network topology and
geometry, hold well for preheated foams. To track how well these laws are
maintained as the foam ages, we introduce metrics for measuring a foam’s
disorder over time and build simple autonomous models for these metrics. While
the quadratic Lewis Law is found to hold well throughout the rupture process,
the Aboav law breaks down rapidly when the Gini coefficient, used for
measuring disparity of cell areas, is approximately .8.

††preprint: APS/123-QED

##  I Introduction

A common topic in materials science is the study of microstructure and its
evolution under various coarsening methods. For planar network microstructure
appearing in polycrystalline metals, porcelains, and foams, the most widely
researched coarsening process occurs through the continuous evolution of grain
boundaries [17, 25, 26, 19, 18]. The coarsening of metals is induced by
annealing, and has a direct relation to its tensile strength [30, 16] . For
foams, coarsening is driven by the transfer of gas between cells with unequal
pressures [40]. In two dimensional planar networks, topological changes are
triggered through two event types. From the Von Neumann-Mullins $n-6$ rule,
cells with fewer than six sides will shrink at constant rates proportional to
their number of sides minus six, eventually shrinking to a point [37, 29].
When this occurs, the network maintains its trivalent structure through
introducing edges in neighboring cells, known as a T2 move. An individual edge
can also shrink to a point and trigger a T1 move which induces topological
changes to its four neighbors [39].

In this study, we study the coarsening process driven exclusively by the
rupture of cell boundaries. This process occurs on a much faster scale than
gas diffusion, and is typically encountered in day to day encounters with
foams. Ruptures produce topological reactions which are markedly different
from those found in T1 and T2 moves. Denoting $C_{n}$ for a cell with $n$
sides, the reaction for the change of topology of the four neighboring cells
$(C_{i},C_{j},C_{k},C_{l})$ bordering a rupturing edge is (typically) given by
the three sub-reactions

$\displaystyle C_{i}+C_{j}\rightharpoonup C_{i+j-4}$ $\displaystyle\textrm{
(Face-merging)},$ (1) $\displaystyle C_{k}\rightharpoonup C_{k-1},\quad
C_{l}\rightharpoonup C_{l-1}$ $\displaystyle\textrm{ (Edge-merging)}.$

The merging of the $i$ and $j$ sided cells is a second order reaction, with
similarities to the sticky particle models of Smoluchowski [33]. Here,
clusters $A_{i}$ and $A_{j}$ with sizes $i$ and $j$ observe the reaction
$A_{i}+A_{j}\rightharpoonup A_{i+j}$. The topological reaction (1) was the
basis for graph and mean-field models of foam rupturing in [20]. In these
models, the reaction rate for face-merging resembled the Smoluchowski equation
with a multiplicative kernel, which is well-known to produce gelation
behavior, or the generation of a massive, infinite-sized cluster[2].

Several previous experimental studies for foams focused on the evolution of
network statistics driven through gas diffusion [8, 34, 35, 6, 31]. The
dynamics of rupturing two-dimensional foams, however, is less studied. An
introductory study was undertaken by Burnett et. al. [4] in which bulbs used
in a light box served as a natural way to heat the foam. For these studies,
the chamber between two plates is filled with a foaming solution and then
flushed with gas. The liquid is allowed to drain to produce a foam of desired
wetness. In our study, which is detailed in Section II, a three-dimensional
foam sample is placed on top of a plate and compressed with another identical
plate. This produces a semi-wet foam which is immediately ready for analysis.
Another heated plate is placed atop the two plates to produce a rupturing
process which typically stabilizes within thirty seconds. This streamlined
process allowed us to obtain multiple samples, a total of 30 experiments with
30 snapshots, taken once per second, producing a dataset of 900 foam
snapshots. With these samples we can now produce standard errors and
confidence bands for the multiple statistics used to analyze foam properties.
To account for open regions near the foam’s border, we impose an artificial
circular boundary which serves as a wall for bordering cells.

In Section III.1, we report on general observations produced from the
experiments. The 30 preheated foam samples have variable initial conditions,
both in total number of cells and distribution of cell areas. This is in
contrast to studies such as [3], which is able to design foams with prescribed
lattice structures through an intricate system of vacuums. We also observe
that rupturing produces massive, irregularly shaped regions, a phenomenon both
observed by Burnett [4] and produced in the computational study of one of the
authors [20]. See Fig. 1 for snapshots of a single foam in the rupturing
process before and after the circular boundary is added. In Section III.2, we
conduct a more detailed analysis of foams before heating. We find uniformity
in certain statistics such as the statistical topology of cells. We also show
strong fits with the Aboav linear law [1] and a quadratic version of the Lewis
law [23], which are empirical laws used for summarizing the topology and
geometry of a planar network.

Figure 1: Snapshots for a foam sample before (images (a), (b), and (c)) and
after (images (d), (e) and (f)) an artificial circular boundary is added.
Left: Initial conditions before heating with 5512 total cells. Middle: The
foam in its intermediate stage after 14 seconds with 3810 total cells. Right:
The aged foam after 29 seconds, with 3222 total cells remaining.

In Section III.3 we model the evolution of a foam’s disorder through
autonomous differential equations. The time scale used in these models is
internal, using the rupture fraction, or the ratio of remaining cell number
over initial cell number. A similar approach is also taken in [4]. For
measuring disparity in cell areas, we use the Gini coefficient [15], a
quantity typically used in measuring income and other economic disparities. We
also present a measure for topological defect through the second moment of
foam’s side number distribution minus six. We find that the Gini coefficient
has more desirable properties at measuring network disorder, and it is against
this measure that we compare how the Aboav and Lewis laws hold as a foam ages.
In Section III.3.3, we find a rapid decay of the correlation coefficient for
the Aboav law when the Gini coefficient is approximately .8. This breakdown
provides us a means for detecting a phase transition in foam disorder, which
in the mean field models of foam rupture in [20] can be more easily found by
simply tracking the side number of the largest cell. The Lewis law, on the
other hand, holds well throughout the rupturing process. In Section IV, we
discuss the implications of these results for modeling foams computationally
and measuring the disorder for foams found in manufacturing and industry.

##  II Experimental methods

A soap solution is created using 1 teaspoon ($\sim 4.93$ cm3) of Palmolive
brand liquid dish detergent and 500 ml of water in a 1 liter container. A lid
is placed on the container and shaken vigorously to create a three dimensional
soap foam which sits atop the liquid solution. A foam sample of approximately
30 cubic centimeters is scooped from the container and placed upon a glass
plate with dimensions of $20\times 25\times.24$ cm3, which is laid on top of a
Cricut LED light box. Spacers of height .5 mm are then added to the corners of
the plate, and a second plate (with same dimensions as the first) is then
slowly placed on top of the first plate. We note that the height of the
spacers are critical to producing foams amenable to analysis. Spaces which are
too tall produce foams with more than a single layer of cells, and those too
short cause issues with overly thin Plateau borders which rupture immediately
when heated. This compresses the three-dimensional foam to form a single-cell,
or quasi two-dimensional, structure which is approximately circular. This foam
sample in this state has a radius of about 15 cm, and there are typically
around 5000 total cells in the network before any heating is applied.

The foam at this point is quite stable, with an occasional rupture occurring
every few seconds. No heat is emitted from the LED light box, so we instead
choose to heat from above by using a heated third plate. This plate is heated
uniformly with a heating pad to $45^{\circ}$C, and is then slowly placed on
top of the two plates enclosing the foam. Recall that spacers are placed
between the plates, so adding the third plate does not further compress the
foam. Nevertheless, the foam does change in size, at first shrinking when
heated, and then slowly expanding as the plate cools. This is due to the
surface tension of the system. For an isotropic soap foam, surface tension is
proportional to the total perimeter of the geometric graph, thereby resulting
in a uniform shrinking of the foam when its temperature rises.

During the heating of the foam, cells shrink but few ruptures occur. Once the
foam reaches its smallest size, cells walls become weak enough from the
heating for the rupture process to begin. The process begins slowly, and then
quickly speeds up to several hundred ruptures per second. For all foams, after
30 seconds has elapsed, the occurrence of ruptures again becomes infrequent.
Snapshots of the process are taken with a Sony Alpha 7 II camera with pixel
resolution $4000\times 6000$ pixels for an image of size $16.8\times 25.2$
cm2. The camera is placed above the foam at a distance of 80 cm. Using an
intervalometer, a snapshot is taken at $t=0,\dots,29$ seconds to give 30
snapshots for the rupturing process, with the first snapshot taken after the
initial shrinking of the foam. This experiment is repeated for 30 different
foam samples, for a total of 900 foam snapshots.

###  II.1 Image processing

In [4], the space between foam bubbles are distinguished from the interior
regions inside bubble walls to define a ‘percentage area covered’. For the
metrics used in Section III.3, large gaps between cells are definitive in
measuring a foam’s disorder, so we will consider all regions as individual
cells. With a free boundary, however, an issue arises when trying to quantify
the size of large, sometimes labyrinthine, regions created from ruptures, but
which are technically part of the foam’s exterior. Our approach for including
these regions is to create an artificial circular boundary. Choosing a circle
as a boundary shape is a natural choice, as the foam maintains a roughly
circular shape throughout the rupturing process. The algorithm of creating the
foam boundary is as follows:

  1. 1.

Compute the centroid $(\bar{x},\bar{y})$ of cell boundary pixels in the image.

  2. 2.

Determine a minimal radius $R$ such that at least 90% of all pixels are
contained in the disk centered at $(\bar{x},\bar{y})$.

  3. 3.

Crop the image by inserting an artificial circular boundary centered at
$(\bar{x},\bar{y})$ with radius $R$ and removing all pixels outside of this
circle.

The regions sharing a border with the boundary are now considered as distinct
cells. Circular arcs on the boundary are counted as cell edges. We have
selected to include 90% of cells since the intersection between the outer
boundary of the foam and the circle is minimal, avoiding the creation of
abnormally large cells which wrap around the boundary.

Several image processing and morphology packages from Python were employed to
process image data from the cropped images. Adaptive thresholding is applied
to the image to account for any imbalances in lighting. The image is then
binarized, followed by an opening (erosion followed by dilation) operation to
remove spurious pixels, and then dilated again to thicken the boundaries
between cells. This prevents against the identification of two cells as a
single connected component. Dilation also prevents the counting of Plateau
border regions found at the foam’s triple junctions. For the semi-wet foams we
consider, the Plateau borders sometimes contain small, Apollonian-like, cells,
which are also coarsened out from dilation. We mention that studies have begun
to consider these ‘inner cells’ in wet foams [14, 32, 22]. Dilating boundaries
increases edge thickness, and subsequently decreases cell areas, but since the
perimeter of the entire foam is uniformly enlarged, this operation has no
discernible effect when comparing relative areas of cells.

Edge detection algorithms then determine connected components of the binarized
image, from which areas can be found readily by a simple pixel count (with a
conversion factor of 1 pixel = $42\times 42\mu$m2). For each cell, the number
of neighbors (or sides) were found by determining connected components (faces)
of the planar network, and constructing an adjacency matrix between the cells.
This enables us to determine first-order topological correlations used in
finding fits for the Aboav law.

##  III Results

###  III.1 General observations and an internal timescale

In Figure 1, we show a sample foam at three stages in its evolution, with and
without the artificially imposed circular boundary. In its initial conditions,
the foam is somewhat uniform in its area and side distributions, although we
still observe several large cells with over fifty neighbors. While gelation
occurs at different times for each foam, we note for this sample foam that by
the middle of the experiment, at 14 seconds, a multitude of large cells have
formed. Merging of both large and small cells continues until the end of the
experiment at 29 seconds, where the interiors of the 20 largest cells comprise
most of the enclosed disk. Over all samples, the most evident observation for
rupturing foams is the rapid generation of multiple large regions. As
mentioned in Section II.1, we consider these regions as proper cells. As time
progresses, these regions grow and are, in general, irregularly shaped and
nonconvex. Surrounding these large regions are smaller, convex cells, with
thickening Plateau boundaries as the liquid content of the foam is distributed
to a smaller amount of cells. These smaller cells often form thin, bridgelike
structures between the massive cells.

The rupture rate for foams depends on its initial wetness. Some foams begin
rupturing almost immediately after their shrinking period, and generate
massive cells within five seconds. Others take several seconds before
beginning a slow rupturing process producing multiple large cells. The timing
for ruptures, in general, has been found to be erratic [36], generally
occurring in cascades with nontrivial spatial correlations. This is in
contrast to coarsening under gas diffusion, where numerical and physical
experiments demonstrate the annihilation of cells at a linear rate which
occurs at uniformly across the foam [9]. A spatial correlation of ruptures in
our experiments also appears to exist, although for foams with a large number
of ruptures, massive cells are approximately uniformly distributed across the
circular domain by the end of the rupturing process.

To use a time scale which is more amenable to dynamic modeling, we will work
with an ‘internal clock’ of the foam, in which we track total cell numbers
relative to initial conditions. In nearly all ruptures, a single rupture
follows the reaction (1), and decreases the total cell count by one. We will
call this time scale the rupture fraction as the total number of ruptures is
approximately the total reduction in cells (see [20] for some counterexamples
where ruptures can remove more than one cell). For a planar network $G_{t}$
and $|G_{t}|$ denoting the total number of cells after $t$ seconds, we define
the rupture fraction as

$s(t)=1-\frac{|G_{t}|}{|G_{0}|}.$ (2)

In Figure 2, we plot $s(t)$ for $t=1,\dots,29$ for each of the 30 foam
samples. Note that in some slower rupturing foams, we observe a ‘time
reversal’ near the beginning where movement of smaller cells near the
artificial circular boundary causes a slight increase in the total number of
cells in the circular region. For the purposes of this study, in which we are
more interested in the behavior for highly aged foams, this effect is minor in
terms of modeling considerations.

Figure 2: The rupture fraction for 30 foam samples, with snapshots taken at
$t=0,\dots,29$ seconds.

###  III.2 Statistics of initial conditions

The three-dimensional foams created from shaking a foam solution vary in
wetness and average cell size. Foam samples are taken from near the
liquid/foam interface, where cell boundaries are thicker and bubble volumes
are smaller. The number of initial cells ranges from 3381 to 10743. A
histogram of total initial cell numbers is given in Fig. 3, and the variation
of wetness produces differing rates of rupture as shown in Fig. 2. While we
find variation among the samples in terms for total initial cell numbers and
wetness, we show in this section we will find that these foams have several
similar network statistics.

Figure 3: A histogram of total cells taken before the heating process for 30
foam samples.

####  III.2.1 Side number and area distributions

The side number distribution, also called the statistical topology, gives
proportions $p_{k}$ of $k$ sided cells for $k\geq 1$. A key observation in
coarsening by diffusion is that under a wide range of ordered and disordered
initial conditions, the statistical topology distribution converges to a
universal attractor, even though the network continues to coarsen with average
cell area increasing at a constant rate [9]. In Fig. 4, we plot the side
distribution for the 30 foam samples before the rupture process begins. As
expected with a cubic planar graph (having all vertices of degree three), the
mean number of neighbors $\langle n\rangle=\sum_{n\geq 1}p_{n}n$ is
approximately 6 , with $\langle n\rangle=5.9\pm.1$ (we report this and all
future confidence intervals with plus or minus one standard error). Despite
the range of total numbers of foams, the side distribution for all of the
foams are similar, with the mode of almost all foams occurring at 4 sides. We
contrast this distribution to the universal attractor distribution from
numerical studies of coarsening driven by mean curvature flow in of [9]. This
distribution (also plotted in Fig. 4) is more concentrated near its mode of 6
sides, and rarely has cells with more than 10 sides.

We also plot the area-weighted side number distribution in Fig. 4, where
$\hat{p}_{k}$ for $k\geq 1$ gives the probability that a randomly selected
interior point in a foam is inside of a $k$-gon. This weighted distribution is
used in finding the average coarsening rate from gas diffusion in [31]. As
cells with more sides tend to be larger (an immediate consequence of the Lewis
law), we should expect, and indeed find, that using area-weighted
distributions increase both the mode and tail probabilities. Denoting the mean
area-weighted number of sides as $\llangle n\rrangle=\sum_{n\geq
1}\hat{p}_{n}n$, we find $\llangle n\rrangle=10.2\pm.9$. Distributions across
different samples are also more erratic than unweighted distributions. In
particular, ‘spikes’ with proportions of several percent frequently arise from
a few large cells all having the same side number.

Figure 4: Side distributions of foams before heating. In both figures,
transparent lines correspond to distributions of individual samples, and the
averages of these distributions are plotted with a bold dashed line.
Connecting lines between integer values of neighbors serve as a visual aid.
Left: Proportion of cells with $k$-neighbors for $k=1,\dots,15$. Also shown
for comparison is the stationary distribution for coarsening through mean
curvature found in [9]. Right: Area-weighted proportion of cells with
$k$-neighbors for $k=1,\dots,30$.

We also plot the distribution of cell areas in Fig. 5. To better visualize the
relation between cell sizes, we scale cells areas to have a mean of 1, and
plot the base ten log of these areas. Even before heating, we find that cell
areas exhibit a multiscale behavior, with large cells with many sides are
surrounded by many smaller 3 and 4-gons with areas differing by several orders
of magnitude. Cell areas appear to vary more than statistical topologies when
considered across different foam samples. However, we find that most
distributions have roughly the same shape—generally left skewed, with a mode
occurring near the mean cell area.

Figure 5: Density of log relative areas of foams before heating. Individual
densities are plotted with transparency, and the mean density over all foams
are shown with a dashed line.

####  III.2.2 Aboav and Lewis Laws

For a cell with $n\geq 2$ sides, Aboav’s law [1] is an empirical observation
that the average number of sides $m_{n}$ for neighboring cells can be
approximated by

$m_{n}=a+\frac{b}{n}\quad\Rightarrow\quad M_{n}:=nm_{n}=an+b.$ (3)

Typically, $M_{n}$ is plotted instead of $m_{n}$, and the model’s fit is then
measured with linear regression. Computing $M_{n}$ is a straightforward
exercise when equipped with the adjacency matrix from the dual graph of the
foam.

With all 30 samples and considering cells with $n=2,\dots,12$ sides, we
compute linear regression parameters for $M_{n}$. We find a strong
correlation, with the average coefficient of determination $r^{2}=.98\pm.01$.
The coefficients of Aboav’s law have values $a=5.4\pm.5$ and $b=23\pm 4$.

Another empirical law related to cellular microstructure is the Lewis law
[23], which proposes a linear relationship between side number $n$ and the
average area $A_{n}$ of $n$-sided cells. For a linear fit $A_{n}=cn+d$ of
relative areas, we find an average coefficient of determination
$r^{2}=.909\pm.003$, but plots for $A_{n}$ appear to be convex, with a
residual analysis suggesting that linear regression is not an appropriate
model. Instead, we consider a quadratic fit
$A_{n}=\beta_{2}n^{2}+\beta_{1}n+\beta_{0}$. See Fig. 6 for $A_{n}$ plotted
over across all samples and the average linear and quadratic Lewis laws. The
coefficient of determination for the quadratic model for each foam sample is
$R^{2}>.999$ with coefficient values $\beta_{0}=.0116\pm.0005$,
$\beta_{1}=-.060\pm.003$, and $\beta_{2}=.080\pm.007$.

Figure 6: Average relative areas $A_{n}$ of $n$-sided cells for $n=2,\dots,12$
over all foam samples plotted with transparency. The linear and quadratic
Lewis law fits use average values of model parameters across all 30 samples.

###  III.3 Dynamic statistics of gelation

As shown in Fig. 1, we find that the second order reaction (1) induced by edge
ruptures creates increasing disparities between massive cells bordering
hundreds of cells and small, mostly unaffected three and four-sided cells. The
reassignments of side number are related to those associated with the
Smoluchowski equation for sticky particle clusters with a multiplicative
reaction kernel. For this model, at a positive finite time a single cluster
experiences an exposive growth producing a massive, infinite sized cluster
called the gel. Borrowing terminology for our study, we will refer to the
creation of massive cells with many sides as gelation. In creating kinetic
models for both cluster models [24, 27] and those for coarsening behavior of
foams [28, 11, 13, 12, 21, 20], all-to-all connectivity is typically assumed.
This assumption leads to gelation behavior which is concentrated in a single
cluster or cell. However, Fig. 1 reveals that gelation in rupturing foams is
often shared between several cells, necessitating metrics which gives similar
results for when a large area is shared by either one or a small number of
massive cells. We develop metrics in this section for measuring gelation from
the perspective of area and topology.

####  III.3.1 Gelation of area

Figure 7: Curve fitting for $\hat{\mathcal{G}}(s)$. Left: Linear regression of
$\hat{\Phi}(g)=ag+b$, with $a=-3.02$ and $b=2.77$. Right: The solution curve
$\hat{\mathcal{G}}(s)$ with initial condition $g_{0}=.6$. For each foam sample
we also overlay, with transparency, the translated gelation curves
$\mathcal{G}_{i}(s-\hat{\mathcal{G}}^{-1}(\mathcal{G}_{i}(0)))$ for
$i=1\dots,30$.

An effective measure for measuring disparities in cell areas is given by the
Gini coefficient, originally formulated to compute the income and other
economic disparities among populations [15]. We will use this measure to
measure disparities among cell areas in a foam. Specifically, let
$A=(a_{1},\dots,a_{N})$ be a nondecreasing sorted list of cell areas. For
$q\in[0,1]$, define $g(q)$ by letting
$g(k/N)=\sum_{i=1}^{k}a_{i}/\sum_{i=1}^{N}a_{i}$ for $k=1,\dots,N$ and
defining other points in $[0,1]$ through linear interpolation. The Gini
coefficient for cell areas is then

$\mathcal{G}=1-2\int_{0}^{1}g(q)dq.$ (4)

For a hexagonal lattice on a fixed circular boundary, $\mathcal{G}$ approaches
zero as the number of cells becomes large. At the other extreme, a single
massive cell encompassing almost the entire domain border by many small cells
can produce a Gini coefficient arbitrarily close to 1.

Using the rupture fraction $s\in[0,1]$ in Eqn. (2) as a time variable, we look
to model gelation as an autonomous differential equation, written as

$\frac{d}{ds}\mathcal{G}(s)=\Phi[\mathcal{G}(s)],\quad\mathcal{G}(0)=g_{0}.$
(5)

In most cases, we find that the Gini coefficient increases in time, so in
general we should require that $\Phi(g)$ is a continuous, positive function
for $g\in[0,g_{*})$ with $g_{*}\leq 1$. As a first model, we will examine the
linear fit

$\Phi(g)=ag+b,\qquad a<0.$ (6)

This produces, for initial conditions $\mathcal{G}(0)=g_{0}$, the solution
curve

$\mathcal{G}(s)=e^{as}\left(g_{0}+\frac{b}{a}\right)-\frac{b}{a}.$ (7)

Under this model, $\Phi(g)$ is positive for values less than $g_{*}=-b/a$, and
this value is an upper bound for the Gini coefficient for all rupturing foams
with initial Gini coefficients $g_{0}<g_{*}$.

A fit for $\hat{\Phi}$ is found with central differences, using our dataset
$(s_{i}^{j},\mathcal{G}_{i}^{j})$ denoting, respectively, the rupture fraction
and Gini coefficient for the $i$th foam taken after $j$ seconds, where
$i=1,\dots,30$ and $j=1,\dots,29$. Derivatives at rupture fractions
$\hat{s}_{i}^{j}=(s_{i}^{j+1}+s_{i}^{j})/2$ are then approximated by

$\mathcal{G}^{\prime}(\hat{s}_{i}^{j})\approx\frac{\mathcal{G}_{i}^{j+1}-\mathcal{G}_{i}^{j}}{s_{i}^{j+1}-s_{i}^{j}}.$
(8)

Fig. 7 shows the scatter plot and a linear fit for the collection of points in
(8). With a linear model, we find a fit of $\hat{\Phi}(g)=-3.02g+2.77$ . From
(7), we use the regression parameters of $\hat{\Phi}$ with $g_{0}=.6$ (chosen
to be smaller than the Gini coefficient for all foams) to obtain the fit
$\hat{\mathcal{G}}(s)$.

Under the time scale of rupture fraction, each of the foams we have observed
start at different initial times. However, we can use the autonomous property
of our model (5) to align times so that the initial conditions of Gini
cofficients all lie on the same solution curve $\hat{\mathcal{G}}(s)$ with
predetermined initial condition $g_{0}$. Specifically, given the initial
condition $\mathcal{G}_{i}(0)=g_{0}^{i}>g_{0}$, we predict

$\hat{\mathcal{G}_{i}}(s)=\hat{\mathcal{G}}(s+\hat{\mathcal{G}}^{-1}(g_{0}^{i})).$
(9)

To visualize how the solution curve $\hat{\mathcal{G}}(s)$ compares against
each of the empirical graphs of Gini coefficient $\mathcal{G}_{i}(s)$, we plot
$\hat{\mathcal{G}}(s)$ against the translated curves
$\mathcal{G}_{i}(s-\hat{\mathcal{G}}^{-1}(g_{0}^{i}))$ so that initial
conditions for empirical curves lie on the solution curve. We find a
reasonable fit between the empirical curves and the model, although for two
foams with smaller initial Gini coefficients we observe a considerable lag in
growth compared to the model.

####  III.3.2 Gelation of topology

We can also measure gelation of a foam from a topological perspective by
computing a distance to an ordered foam containing only hexagons. To define
this measure, we consider the random variable $S$ of the side number of a
random cell selected with uniform probability. For the probability mass
function $\\{p_{s}\\}_{s\geq 1}$ denoting probabilities of selecting $s$-sided
cells, we define the topological gelation which measures the second moment of
the topological defect $S-6$, given by

$\mathcal{T}=\mathbb{E}[(S-6)^{2}]=\sum_{s\geq 1}p_{s}(s-6)^{2}.$ (10)

Like the Gini coefficient, $\mathcal{T}$ is nearly 0 for a perfect hexagonal
lattice on a circular boundary. However, unlike the Gini coefficient,
topological gelation can be arbitrarily large. In particular, a foam
consisting of a single interior cell of $n$ sides adjacent to a band of $n$
4-sided boundary cells gives a topological defect that is asymptotic to $n$ as
$n\rightarrow\infty$.

Topological gelation appears to grow faster than a linear rate, and so to
check for power law growth we compute linear fits for log data
$\widehat{\ell}_{i}:=\log_{10}(\widehat{\mathcal{T}}_{i}(s))=as+b$ for each
observed gelation curve $\mathcal{T}_{i}(s)$. The average of these
coefficients is computed to give a fit $\widehat{\ell}(s))=as+b$ with $a=2.22$
and $b=1.14$. We plot this regression line in Fig. 8. Similar to shifting
gelation curves for the Gini coefficient, we assume the model is autonomous
and plot $\ell_{i}(s-\widehat{\ell}^{-1}(\ell_{i}(0)))$ against the model
$\widehat{\ell}(s)$. It appears, however, that the growth rate in topological
defect increases more highly aged foams. For foams with $\mathcal{T}<100$,
growth scales at approximately $O(s^{2})$, while for $\mathcal{T}>100$, we
find a growth closer to $O(s^{4})$. We also note that a few highly aged foams
sometimes decrease in $\mathcal{T}$. This can be an effect of several
artifacts, including the movement of foam pockets into the artificial
boundary, and the fusing of two regions of bubbles to reduce total variance of
cell neighbors.

Figure 8: The linear fit $\widehat{\ell}(s)=as+b$ with $a=2.22$ and $b=1.14$
for the log of topological gelation. For each foam sample we also overlay,
with transparency, the translated gelation curves
$\ell_{i}(s-\widehat{\ell}^{-1}(\ell_{i}(0)))$ for $i=1\dots,30$.

####  III.3.3 Decay of fit for Aboav and Lewis laws

In Section III.2.2, for unheated foams we found a strong coefficient of
determination value for both the linear Aboav law and quadratic Lewis law. As
the foam ages, it is not clear if the parameters in these models remain
constant, or even if these models keep their high values of $R^{2}$. From
Section III.3, we find that the growth of the Gini coefficient is smoother and
more predictable than topological gelation, so this will be the metric we will
compare against the Aboav Law. In Fig. 9, we plot $\mathcal{G}$ against the
correlation coefficient over all snapshots. We also plot $\mathcal{G}$ against
the fitted parameter values $a$ and $b$ in the Aboav fit $M_{n}=an+b$. To
create smoothed estimates we use locally-weighted scatterplot smoothing
(LOWESS) regression [7]. To create a 95% confidence band, we generate 100
bootstrap samples of LOWESS curves, each derived from randomly selecting half
of the data.

Figure 9: Scatterplots and LOWESS regression curves for Aboav’s Law
$M_{n}=an+b$ over all foam snapshots, comparing the Gini coefficient
$\mathcal{G}$ with (a) correlation coefficient, (b) value for intercept
parameter $b$, and (c) value for slope parameter $a$. For each image, the 95%
bootstrap confidence band is shown in light grey.

We find that the linear Aboav fit holds quite well when $\mathcal{G}<.8$, with
$r=.975\pm.001$. Fitted parameter values are also stable in the regime, with
$a=5.28\pm.03$ and $b=28.4\pm.5$. When $\mathcal{G}\geq.8$ the Aboav fit
begins to quickly decay, with the average correlation reducing to
$r=.35\pm.03$. The slope parameter $a=4.9\pm.4$ decreases somewhat but becomes
more variable as $\mathcal{G}$ increases, while the intercept $b=270\pm 20$
grows rapidly, due to the high probability of low-sided grains neighboring
massive cells.

Recall that the Lewis law does not consider area or side correlations for
neighboring cells, but rather simply provides average relative areas of
$n$-gons. The generation of massive cells appears to have little effect on the
average relative areas on neighboring smaller cells. This is observed in
Figure 10, where we consider relative area data in aggregate for each
$n$-sided cell for $n=2,\dots,12$. Variance of relative areas grows between
samples as the number of sides increases, with the largest standard error of
.03 occurring for 12 sided cells. For each $n=2,\dots,12$, we aggregate
relative areas across all images, and then compute mean relative areas. The
resulting graph fits well against the quadratic Lewis law
$A_{n}=\beta_{2}n^{2}+\beta_{1}n+\beta_{0}$, with coefficient of determination
$R^{2}=.99$ and fitted values $\beta_{2}=.011$, $\beta_{1}=-.063$, and
$\beta_{0}=.094$.

Figure 10: Average relative areas $A_{n}$ of $n$-sided cells for
$n=2,\dots,12$ over all foam samples, plotted with standard errors. The
quadratic Lewis law fit $A_{n}=.011n^{2}-.063n+$ using all data points is
shown in the dashed line.

##  IV Discussion

In this study, we have gathered statistics for a class of two-dimensional
rupturing soap foams. The method for generating and rupturing foams is simple
to reproduce. However, given the stochastic nature of the rupturing process we
find variation in area and side number statistics among both heated and
unheated foams. Because we were able to replicate the experiment for 30
different foam samples, we provided error estimates for our statistics. The
foams in this study were restricted to soap foams, but the gelation measures
introduced in this work could be applied to any dynamic two-dimensional planar
network. It would be interesting to investigate if experimental parameters
such as viscosity could affect spatial correlations of ruptures and,
subsequently, statistics such as the Aboav law.

Like the initial studies of the Aboav and Lewis laws, the findings in this
paper are empirical, but theoretical models and analytic tools might help
explain why they are good fits. In particular, a deeper investigation into
gelation behavior might elucidate why a breakdown of the Aboav Law appears
when the Gini coefficient becomes large, and whether this breakdown occurs
simultaneously with a gelation time in a kinitic model. In future work, we
hope to compare experimental data on side-distributions against the
computational rupture model in [20] which focused only on topological
gelation. To incorporate areas in kinetic models of foam evolution, we could
generalize equation (1) to include the merging of cell areas. However, there
exists no known analog of the $n-6$ rule for area evolution in rupturing
foams, as the dynamics of comparing cells boundaries before and after rupture
is complicated. One simple approximation may consider an additive model for
areas, in which a rupture causes two cells of areas $A_{1}$ and $A_{2}$ to
merge into a single cell of area $A_{1}+A_{2}$, with all other cells
maintaining their original areas. For a foam comprised of cell sides and areas
$\\{(s_{i},a_{i})\\}_{i=1}^{N}$, the reaction resulting from a rupture is then

$\displaystyle(s_{i},a_{i})+$
$\displaystyle(s_{j},a_{j})+(s_{k},a_{k})+(s_{l},a_{l})$ (11)
$\displaystyle\rightharpoonup(s_{i}+s_{j}-4,a_{i}+a_{j})+(s_{k}-1,a_{k})+(s_{l}-1,a_{l}).$

The limiting integro-differential equation for the evolution of area and side
distributions would be unusual in that it would combine a continuous merging
for areas with discrete merging for side-numbers. A discrete time stochastic
particle system following (11) can compare statistics against those found in
this study. The Gini coefficient is particularly amenable to such a reduced
model, as it only requires a sorted list of cell areas to compute.

Assuming no memory effects, the evolution of cell areas and topologies can be
modeled as a homogeneous Markov process. For a state space consisting of a
list of cell sides and areas, transitions would be determined by randomly
selecting four cells to undergo the reaction (11). In most mean field models
for grain coarsening, cells are chosen solely in proportion to their number of
sides [28, 11, 12] (for a notable exception, Chae and Tabor [5] chose cells
proportional to cell length). In [21], for instance, it was shown that adding
first-order topological neighbor correlations in grain coarsening models
produced significant differences in statistical topologies. Incorporating
other topological and geometrical relations into transition probabilities may
help produce more accurate models of network statistics. The Aboav laws, which
in Section III.3.3 were found to hold reasonably well for foams with Gini
coefficient less than .8, can perhaps be used, although other models will
likely be necessary after entering into the non-Aboav, or gelation, regime.

We close by mentioning a class of related coagulation equations which may be
amenable to rigorous analysis. The reaction (1) can be seen as an instance of
a reaction with dissipation, since six sides are lost with each rupture. In
[38], merging reactions for for the Smoluchowski equation were paired with an
additional annihilation reaction $C_{j}\rightharpoonup\emptyset$ for clusters
size $j\geq 1$, and several exact and approximate formulas were found for
distributions of $C_{n}$. An analytic study of merging clusters

$C_{i}+C_{j}\rightharpoonup C_{i+j-\phi(i,j)},\qquad i,j\geq 1$ (12)

for some reaction cost $\phi(i,j)$ remains to be seen. As in [38], a key
question is whether the system completely dissipates before producing a gel.

###### Acknowledgements.

The work of JK and EO is partially supported by a National Aeronautics and
Space Administration, Pennsylvania Space Grant Consortium grant 80NSSC20M0097,
Subaward S002967. The work of JK is partially supported by a National Science
Foundation LEAPS-MPS grant 2316289.

## References

  * [1] D. Aboav, The arrangement of grains in a polycrystal, Metallography, 3 (1970), pp. 383–390.
  * [2] D. J. Aldous, Deterministic and stochastic models for coalescence (aggregation and coagulation): a review of the mean-field theory for probabilists, Bernoulli, (1999), pp. 3–48.
  * [3] J. Bae, K. Lee, S. Seo, J. G. Park, Q. Zhou, and T. Kim, Controlled open-cell two-dimensional liquid foam generation for micro-and nanoscale patterning of materials, Nature communications, 10 (2019), pp. 1–9.
  * [4] G. Burnett, J. Chae, W. Tam, R. M. De Almeida, and M. Tabor, Structure and dynamics of breaking foams, Physical Review E, 51 (1995), p. 5788.
  * [5] J. Chae and M. Tabor, Dynamics of foams with and without wall rupture, Physical Review E, 55 (1997), p. 598.
  * [6] A. T. Chieco and D. J. Durian, Experimentally testing a generalized coarsening model for individual bubbles in quasi-two-dimensional wet foams, Physical Review E, 103 (2021), p. 012610.
  * [7] W. S. Cleveland, Robust locally weighted regression and smoothing scatterplots, Journal of the American statistical association, 74 (1979), pp. 829–836.
  * [8] J. Duplat, B. Bossa, and E. Villermaux, On two-dimensional foam ageing, Journal of fluid mechanics, 673 (2011), pp. 147–179.
  * [9] M. Elsey, S. Esedoglu, and P. Smereka, Large-scale simulation of normal grain growth via diffusion-generated motion, Proceedings of the Royal Society A: Mathematical, Physical and Engineering Science, 467 (2011), pp. 381–401.
  * [10] P. J. Flory, Molecular size distribution in three dimensional polymers. i. gelation1, Journal of the American Chemical Society, 63 (1941), pp. 3083–3090.
  * [11] H. Flyvbjerg, Model for coarsening froths and foams, Physical Review E, 47 (1993), p. 4037.
  * [12] V. Fradkov, A theoretical investigation of two-dimensional grain growth in the ‘gas’ approximation, Philosophical Magazine Letters, 58 (1988), pp. 271–275.
  * [13] V. Fradkov, D. Udler, and R. Kris, Computer simulation of two-dimensional normal grain growth (the ‘gas’ approximation), Philosophical Magazine Letters, 58 (1988), pp. 277–283.
  * [14] N. Galvani, M. Pasquet, A. Mukherjee, A. Requier, S. Cohen-Addad, O. Pitois, R. Höhler, E. Rio, A. Salonen, D. J. Durian, et al., Hierarchical bubble size distributions in coarsening wet liquid foams, Proceedings of the National Academy of Sciences, 120 (2023), p. e2306551120.
  * [15] C. Gini, On the measure of concentration with special reference to income and statistics, colorado college publication, General series, 208 (1936).
  * [16] E. Hall, The deformation and ageing of mild steel: Iii discussion of results, Proceedings of the Physical Society. Section B, 64 (1951), p. 747.
  * [17] T. Ilmanen, A. Neves, F. Schulze, et al., On short time existence for the planar network flow, Journal of Differential Geometry, 111 (2019), pp. 39–89.
  * [18] L. Kim and Y. Tonegawa, On the mean curvature flow of grain boundaries, in Annales de l’Institut Fourier, vol. 67, 2017, pp. 43–142.
  * [19] D. Kinderlehrer and C. Liu, Evolution of grain boundaries, Mathematical Models and Methods in Applied Sciences, 11 (2001), pp. 713–729.
  * [20] J. Klobusicky, Markov models of coarsening in two-dimensional foams with edge rupture, Journal of Nonlinear Science, 31 (2021), pp. 1–27.
  * [21] J. Klobusicky, G. Menon, and R. L. Pego, Two-dimensional grain boundary networks: stochastic particle models and kinetic limits, Archive for Rational Mechanics and Analysis, (2020), pp. 1–55.
  * [22] S. Kwok, R. Botet, L. Sharpnack, and B. Cabane, Apollonian packing in polydisperse emulsions, Soft Matter, 16 (2020), pp. 2426–2430.
  * [23] F. T. Lewis, The effect of cell division on the shape and size of hexagonal cells, The Anatomical Record, 33 (1926), pp. 331–355.
  * [24] A. A. Lushnikov, Some new aspects of coagulation theory, Izv. Akad. Nauk SSSR, Ser. Fiz. Atmosfer. I Okeana, 14 (1978), pp. 738–743.
  * [25] A. Magni, C. Mantegazza, and M. Novaga, Motion by curvature of planar networks ii, arXiv preprint arXiv:1301.3352, (2013).
  * [26] C. Mantegazza, M. Novaga, and V. M. Tortorelli, Motion by curvature of planar networks, Annali della Scuola Normale Superiore di Pisa-Classe di Scienze, 3 (2004), pp. 235–324.
  * [27] A. H. Marcus, Stochastic coalescence, Technometrics, 10 (1968), pp. 133–143.
  * [28] M. Marder, Soap-bubble growth, Physical Review A, 36 (1987), p. 438.
  * [29] W. W. Mullins, Two-dimensional motion of idealized grain boundaries, Journal of Applied Physics, 27 (1956), pp. 900–904.
  * [30] N. J. Petch, The cleavage strength of polycrystals, J. Iron Steel Inst., 174 (1953), pp. 25–28.
  * [31] A. Roth, C. Jones, and D. J. Durian, Bubble statistics and coarsening dynamics for quasi-two-dimensional foams with increasing liquid content, Physical Review E, 87 (2013), p. 042304.
  * [32] S. Sauerbrei, E. Haß, P. Plath, et al., The apollonian decay of beer foam bubble size distribution and the lattices of young diagrams and their correlated mixing functions, Discrete Dynamics in Nature and Society, 2006 (2006).
  * [33] M. v. Smoluchowski, Drei vortrage uber diffusion, brownsche bewegung und koagulation von kolloidteilchen, ZPhy, 17 (1916), pp. 557–585.
  * [34] J. Stavans, Temporal evolution of two-dimensional drained soap froths, Physical Review A, 42 (1990), p. 5049.
  * [35] , Evolution of two-dimensional cellular structures: The soap froth, Physica A: Statistical Mechanics and its Applications, 194 (1993), pp. 307–314.
  * [36] N. Vandewalle and J. Lentz, Cascades of popping bubbles along air/foam interfaces, Physical Review E, 64 (2001), p. 021507.
  * [37] J. Von Neumann, Discussion: grain shapes and other metallurgical applications of topology, Metal Interfaces, (1952).
  * [38] J. A. Wattis, D. G. McCartney, and T. Gudmundsson, Coagulation equations with mass loss, Journal of engineering mathematics, 49 (2004), pp. 113–131.
  * [39] D. Weaire and N. Rivier, Soap, cells and statistics—random patterns in two dimensions, Contemporary Physics, 25 (1984), pp. 59–99.
  * [40] D. L. Weaire and S. Hutzler, The physics of foams, Oxford University Press, 2001.