Stability maps for columnar structures

ABSTRACT

We have previously explored the hysteresis and reversibility of transitions between ordered packings of soft spheres of diameter d in cylindrical channels of diameter D [A. Mughal, J. Winkelmann, D. Weaire, S. Hutzler, Phys. Rev. E 98, 043303 (2018)]. Here we extend these initial results to include transitions between all columnar structures without inner spheres (i.e. packings in which all of the spheres are in contact with the cylindrical boundary). These results can be represented by a directed network showing permissible transitions between structures. From the hard sphere limit we deduce that there are two different types of transitions, reversible and irreversible. We explore the nature of these transitions for soft spheres as a function of pressure and due to changes in the ratio D/d. These results are illustrated by the use of schematic diagrams, indicating the topological features of each transition. Specific cases are tabulated and can be understood by reference to the appropriate schematic diagram.

KEYWORDS:

• Columnar structures
• sphere packing
• hysteresis
• soft-matter

1. Introduction

Ordered dense packings of spheres in cylindrical confinement have attracted increasing attention in recent years, with examples being found at various length scales [1]. On the microscopic scale, many examples of biological microstructures (e.g. microtubules, viruses and flagella) can be described as tubular packings of spheres [2, 3], while on the macroscopic level the arrangement of buds and leaves in plants can be understood by similar means [4].

Columnar structures have also been realised in a range of experiments, some of the earliest include the packing of bubbles inside tubes [5–12]. More recently, it has been shown that colloidal particles also crystallise into microscopic helical structures when trapped inside cylindrical pores [13–16]. In nano-science, similar structures have been generated by inserting molecules or atoms into cylinders [17–20] or on their surface [21] in order to create novel materials. Further elaborations on the same theme include, the ground state arrangement of repulsive point-like particles confined to the surface of a cylinder [22], the packing of spheroids in cylinders [23] and the packing of ellipses on the surface of a cylinder [24]

We and others have published extensive results for the packing of hard [3, 25–28] and soft [12, 16, 29] spheres. The latter are modelled by purely repulsive harmonic forces, allowing for an overlap between spheres. This is a rough representation of bubbles in particular [30]. These results are contained in an extensive phase diagram (Figure 2 in [29]) which specifies the equilibrium structures for all relevant values of D/d (i.e. the ratio of cylinder diameter D to sphere diameter d) and applied pressure p.

Figure 2. Computational results for the densest columnar sphere packings. Shown is the packing density Φ as a function of the ratio D/d of tube to sphere diameter, together with the phyllotactic notation of the structures (see [27] for an extended range of D/d values). The shaded region of the diagram is of relevance to the soft sphere packings discussed in this article. Whenever the line-slip packings $(3,\mathbf{2},1)$ and $(3,\mathbf{3},0)$ are not those of highest density their packing fraction is indicated by dotted lines (these are for illustrative purposes only and not the result of computation).

For the macroscopic systems for which this work has proved relevant so far (foams [6, 7, 31, 32]), questions of metastability and hysteresis arise. The system must remain in a (meta)stable state, but not necessarily that of lowest enthalpy. We have broached this question in our earlier investigation, but provided only one limited example of what we term a ‘stability diagram’, indicating the structures that are encountered for a given change of pressure p or size ratio D/d, and the loci for instability and structural transitions.

In the present paper we will provide a more comprehensive answer to the question raised, by investigating all parts of the phase diagram (up to the point at which internal spheres are found in the structures of maximal density, i.e. $D/d=2.713$ [27]). There are several such cases: we condense the results by representing various types of stability diagrams by a few key parameters.

The paper is organised as follows. In Section 2 we summarise some of our previous results for hard and soft sphere packings. In Section 3 we present the key features of a stability diagram. In Section 4 we introduce a new directed network graph of possible transitions between packings. In Section 5 we extend our results to a number of reversible transitions and also present stability diagrams for irreversible transitions. We give our conclusions in Section 6.

2. Summary of previous results

2.1. Hard sphere packings

For close packings of hard spheres in cylindrical tubes, such as those shown in Figure 1, the key parameter is the ratio $D/d$ between the tube diameter D and sphere diameter d [26, 27]. Varying D/d gives rise to a family of helical structures for which we have provided a complete taxonomy (in terms of phyllotactic indices, see [26, 27] for details) up to the above mentioned critical value $D/d=2.713$. Figure 2 shows the variation of packing density Φ (volume fraction of spheres) with D/d for most of the structures discussed in this paper.

Figure 1. An example of a $(4,4,0)$ uniform (a) columnar structure and its line slip modification (b), generated by packing spheres of diameter d into a tube of diameter D. In the uniform structure, every sphere has an identical neighbourhood with a coordination number z = 6. Its corresponding line slip structure (shown to the right) is characterised by gaps (i.e. loss of contacts), highlighted by the (red) ellipse. For the structure classification we use the phyllotactic notation (see [26, 27] for details).

For a set of discrete values of D/d the densest packing arrangements are given by uniform structures, that is packings in which each sphere is in contact with six neighbouring spheres. For all other values of D/d the densest packings are line-slip arrangements whereby neighbouring spiral chains slide over each other, so that some spheres only have five contacts. An example of a uniform and its related line-slip structure is given in Figure 1.

In Figure 2 the uniform structures correspond to the local peaks in the packing density. Uniform packings can be continuously transformed into the appropriate line-slip arrangement by the loss of a contact, the result is a structure with a lower density. Local minima, on the other hand, are where a given line-slip is surpassed in density by another line-slip structure. The two structures are unrelated to each other, that is, there is no simple continuous transformation between them.

In what follows we will have reason to refer to the region highlighted in grey in Figure 2. For now we simply note that the diagram presents the packing arrangement (either a uniform arrangement of a line-slip packing) with the greatest density for any given D/d, alternative packings with lower density are not shown. However, in the grey region of Figure 2 we illustrate the continuation of some of these line-slip structures (of a lower density) by dotted lines. These dotted lines have been simply sketched in for illustrative purposes and are not the result of computation.

2.2. Soft sphere packings: phase diagram

In [12, 29] we have extended our study of columnar sphere packings to packings of soft spheres. In this case the pressure P is an additional parameter. Stable structures are found by minimising the enthalpy H = E + PV for a system of N soft spheres contained in a volume V and given diameter ratio D/d; E is the internal energy due to sphere-sphere and sphere-cylinder overlaps. The results were presented in the form of an equilibrium phase diagram, showing the structure of minimal enthalpy for a given value of D/d and P [12, 29] (see Appendix A for details of how the simulations were conducted). Boundaries between phases can be continuous or discontinuous, indicating a smooth or abrupt change in the enthalpy, respectively.

For illustrative purposes we present a part of this phase diagram as a schematic in Figure 3 (see Figure 2 in [29] for the full version). The various regions correspond to different lowest enthalpy structures. Their labelling, U or L concern their classification as uniform or line-slip structures, respectively. The regions are separated by phase boundaries at which the change of structure is either continuous (dashed lines) or discontinuous (solid lines).

Figure 3. A schematic version of a part of the full phase diagram presented in [29], showing the structure of minimal enthalpy structure. The horizontal axis (i.e. the pressure p = 0 (hard-sphere) limit) corresponds to the shaded part of Figure 2. The red diamond indicates the value of D/d at which the uniform packing $(4,2,2)$ is found in the hard sphere limit. The regions of the phase diagram indicates which packing – either uniform structures ($U_0,U_1,U_2$) or line-slip ($L_0,L_1,L_2$) has the lowest enthalpy.

The nature of these boundaries can be deduced using the fact that the hard-sphere results (for the shaded region) in Figure 2 correspond to the horizontal axis (P = 0) in Figure 3.

To illustrate the relationship between hard sphere results and the phase diagram for soft sphere packings we shall concentrate on the grey (shaded) region shown in Figure 2. Note in particular that the uniform packing $(3,3,0)$ can be transformed into the line-slip $(3,\mathbf{3},0)$ and this eventual becomes the uniform packing $(4,3,1)$. That is, in the hard sphere limit, as shown by the dotted line, the line-slip $(3,\mathbf{3},0)$ connects the uniform packings $(3,3,0)$ and $(4,3,1)$ by a continuous transformation. Even if for part of the shaded region shown in Figure 2 the line-slip $(3,\mathbf{3},0)$ has a lower density as compared to other packings, nevertheless it still provides a metastable solution. This is reflected in the soft-sphere results, as sketched in Figure 3, whereby the line-slip $(3,\mathbf{3},0)$ continues as a metastable solution (of higher enthalpy) between the two regions labelled $L_1$.

From the hard sphere results it can be seen that the uniform packing $(3,3,0)$ can be transformed continuously into the line slip structure $(3,\mathbf{3},0)$ by the loss of a contact. This relationship extends into the soft-sphere results, as indicated by the dashed line between $U_1$ and $L_1$. As expected, the enthalpy varies smoothly and continuously across the boundary [12].

Similarly, it can be seen from hard sphere results that the uniform structure $(4,3,1)$ can be transformed into $(3,\mathbf{3},0)$ by the loss of a contact, this is again reflected in the phase diagram, Figure 3. In a similar manner, the correspondence between $U_0$ and its associated line-slips ($L_2$ & $L_3$) can also be easily deduced by comparing Figures 2 and 3.

At higher pressures the line-slip arrangements become unstable, vanishing at triple-points, leaving only uniform packings. As such the boundaries between these packings (at high pressures) are given by solid lines to indicate that the change is discontinuous. That is, the change in structure is abrupt and this can also be seen in the enthalpy, which upon crossing the boundary may involve a sudden change in the gradient and (or) a jump-discontinuity (see for example Figure 3 in [12]).

Finally, there are the boundaries between adjacent line-slip packings in Figure 3 – such as the boundary between $L_1$ and $L_2$, or between $L_1$ and $L_3$. These are discontinuous (as expected) since there exists no simple continuous transformations leading from one to the other.

3. Metastability and hysteresis

Given a phase diagram of the type shown in Figure 3, one may ask the question what happens to a given equilibrium structure if the pressure or the diameter ratio is continuously varied? (Experimentally, the first such variation might be easier to realise.) Could the system remain in a metastable state, giving rise to hysteretic effects? Assuming the pressure is sufficiently low, is it possible to begin with the structure $U_1$ and gradually increase the diameter of the confining cylinder so that it transforms first into an intervening line-slip arrangement and then into $U_2$? Could this process be reversed?

The answer is yes, and we have demonstrated this in [29] for the case $U_1=(3,2,1)$ which can be transformed into $U_2=(4,2,2)$ (and back again). We summarised our results in the form of a stability diagram, which we will discuss in Section 5, as it is of relevance to all reversible transitions encountered for the range of D/d analysed here.

Below a critical pressure $p<p_3$ the transformations between uniform and line-slip structures is continuous and is accompanied by a smooth variation in the enthalpy, or at most a change in the derivative of the enthalpy when a new contact is formed (see Appendix A for details of how the simulations were conducted).

The situation is more complex at higher pressures where transitions between structures become discontinuous, leading to hysteresis. Nevertheless, a transformation between the two uniform structures $U_1$ and $U_2$ is possible.

A change in D/d can potentially transform a uniform structure into a number of intermediate line-slip arrangements, each leading to a different uniform structure. In our numerical results we do not observe all of these potential transitions, only a subset, which we call favourable transitions. For a given uniform structure the favourable transition can easily be identified from our previous hard sphere results: it corresponds to the line-slip immediately adjacent to a given uniform structure; see Table 1 in [27] for a list.

Table 1. Critical pressures for reversible transitions.

Structural transition $U_1\iff U_2$	$p_1/{10}^2$	$p_2/{10}^2$	$p_3/{10}^2$
$(2,2,0)\iff (3,2,1)$	0.79	0.40	0.21
$(3,2,1)\iff (4,2,2)$	2.74	2.48	1.47
$(4,3,1)\iff (5,3,2)$	4.70	2.34	0.90
$(3,3,0)\iff (4,3,1)$	3.54	0.87	0.24
$(5,3,2)\iff (6,3,3)$	4.36	2.94	1.18
$(5,4,1)\iff (6,4,2)$	5.36	2.00	0.12
$(4,4,0)\iff (5,4,1)$	5.46	1.02	0.08
$(5,5,0)\iff (6,5,1)$	6.18	0.94	0.08

Note: The values for the pressures $p_1$, $p_2$ and $p_3$ determine the type of transitions encountered upon a change of the size ratio D/d in the stability diagram of Figure 5. This characterises all reversible cases for structures without inner spheres ($D/d<2.713$). The uniform structures $(6,3,3)$ and $(6,4,1)$, which appear above the limit $D/d<2.713$, have also no inner spheres in our cases. All stability diagrams were obtained with a resolution in pressure of $\mathrm{\Delta }p=4\cdot {10}^{-4}$, which determines the uncertainties for all pressures.

This table lists the densest sphere packing structures in sequence with increasing D/d, thus for example the uniform structure $(3,3,0)$ is listed as the 10th structure. Immediately adjacent to this in the table are two line-slip structures: the 9th structure $(\mathbf{3},2,1)$ and the 11th structure $(3,\mathbf{3},0)$, which occur at slightly lower or higher values of D/d respectively. Starting with the uniform structure $(3,3,0)$ and increasing D/d implies that it will transform into the line-slip $(3,\mathbf{3},0)$ which eventually leads to the uniform structure $(4,3,1)$, while starting with $(3,3,0)$ and decreasing D/d implies that it will transform into the line-slip $(3,\mathbf{2},1)$ which eventually leads to the uniform structure $(3,2,1)$.

These are the favourable transitions leading away from a uniform structure, that is, the line-slip with the greatest density the uniform structure can transform into due to an infinitesimal change in D/d.

The favourable transitions (i.e. the line slips) can also be identified from Figure 2. They occur via the line slips immediately adjacent to a local peak (i.e uniform packing) corresponding to the line-slips with the greatest (local) density. Thus for example, starting with the uniform structure $(3,3,0)$ an increase in D/d will transform it into the line slip $(3,\mathbf{3},0)$ and eventually into the uniform structure $(4,3,1)$ – that is we follow the dotted line in Figure 2 that leads from $(3,3,0)$ to $(4,3,1)$ via $(3,\mathbf{3},0)$. If the pressure is below $p_3$ this transformation is continuous, while above this pressure there is hysteresis.

4. Network of structural transitions in the hard sphere limit

To summarise the structural transitions that are encountered upon a change of the ratio D/d of tube to sphere diameter we have developed a directed network graph which indicates a pattern for favourable transitions amongst all possible transitions. From this network graph, shown in Figure 4, we identify two different transition types (reversible and irreversible), as described below. In addition to providing a new perspective on the previous hard sphere results, the network also predicts structural transitions yet to be realised experimentally.

Figure 4. Directed network displaying possible structural transitions between uniform structures, in the hard sphere limit. The structures are restricted to packings without internal spheres and are labelled using phyllotactic notation. The point $(0,0,0)$ represents a conceptual origin for the diagram, corresponding to $D/d=0$. Dashed lines are contours of constant D/d – a transition from one structure to another which involves moving along the diagram in a direction perpendicular to the contours implies a change in D/d. Favourable transitions are marked with blue arrows, unfavourable transitions with red arrows and all unstable cases are shown with black dashed arrows.

The network has as its nodes all uniform structures (specified by their respective phyllotactic notation), arranged so that their distance to some conceptual origin $(0,0,0)$ increases with their values of D/d. The different types of transitions, as detailed below, are indicated by different types of arrows. A transition from one uniform structure to another which requires moving along the diagram in the direction perpendicular to the contours of constant D/d (shown as dashed lines) implies a change in D/d.

Favourable transitions are indicated by bold blue arrows in Figure 4. Also shown, by the red arrows, are the transitions identified from previous hard sphere results [27] which are found to be unfavourable (i.e. intermediate line slip structures exist between these uniform packings but we do not observe these transitions in our numerical simulations). The black dashed arrows indicate line slip transitions that are mechanically unstable. These are of two types as described below.

The first type are rhombic arrangements. That is the transition between packings involves an affine shear of the entire packing (as opposed to a line-slip where the shear is restricted to single ‘line’) which we showed previously to be unstable [33]. In Figure 4 these are transitions between adjacent packings with indices $(n+2,n+1,1)$, i.e. $(2,1,1)\iff (3,2,1)\iff (4,3,1)\iff (5,4,1)\iff (6,5,1)$.

The other mechanically unstable transitions are due to so called re-entrant line slips [26, 27, 33]. These are line-slips for which the cylinder diameter does not increase (or decrease) monotonically with the line-slip. See for example the line slip between $(4,2,2)$ and $(4,3,1)$ in Figure 10 of [27], where the line slip appears to ‘double back’.

Reversible transitions take place between structures connected via two blue arrows, for instance $(3,2,1)\iff (4,2,2)$. When increasing D/d starting from $(3,2,1)$, the preferred transition (blue arrow) is into $(4,2,2)$; similarly upon decreasing D/d the favoured transition returns the system back to $(3,2,1)$ structure.

An example of an irreversible transition is that from $(4,2,2)\Rightarrow (5,3,2)$. Here, the $(5,3,2)$ structure can not be transformed back into the $(4,2,2)$ via decreasing D/d, instead it transforms (reversibly) into $(4,3,1)$.

Following the favourable transitions in the network, we can identify cyclic patterns that can only be followed counterclockwise, with the first cycle starting in the third row: $(3,2,1)\iff (4,2,2)\Rightarrow (5,3,2)\iff (4,3,1)\iff (3,3,0)\Rightarrow (3,2,1)$. The transitions are achieved via increasing D/d until the $(5,3,2)$ structure is reached, at which point decreasing D/d leads back to the original structure via $(4,3,1)$ and $(3,3,0)$. Another potential cycle could start and end with the structures $(5,3,2)$ (confirmation would require further simulations that we have not conducted here).

The structural transitions presented in the network of Figure 4 is broadly consistent with previous results for similar systems. The earliest detailed experimental work appears to be that by [7] who dilated and compressed crystalline arrangements of equal volume soap bubbles that were confined in a cylinder, by the use of a piston. The observed cycles of structural transitions were interpreted by the authors with reference to the nucleation and motion of a dislocation on a triangular lattice [6, 7]. More recently, transitions between assemblies of hard particles in cylinders were explored numerically using Monte Carlo simulations [14]. By compressing the system it is possible to induce transitions between structures, again we find that the favourable transitions listed by Fu et al. [14] are consistent with Figure 4.

5. Stability diagrams for columnar structure with $D/d<2.713$

Having defined the two different transition types (reversible and irreversible) in the previous section, we now want to explore how these transitions change with increasing pressure p.

We do this by computing trajectories in the $(p,D/d)$ plane, and recording boundaries where a structure changes to one of a different character. With a sufficient number of trajectories, a stability diagram is built up, that is, a map of the location of structural transitions.

We will present schematic stability diagrams for both reversible and irreversible transitions, applicable to all the transitions between columnar structures with size ratio $D/d<2.713$ (without inner spheres). The corresponding critical pressures, which demarcate different transition routes, are presented in a table.

5.1. The reversible case

In [29] we had presented, as an example of a reversible transition, the case $(3,2,1)\iff (4,2,2)$. The corresponding stability diagram can be represented schematically, as in Figure 5, so that it shows the topological features of any reversible transition between two uniform structures $U_1$ and $U_2$.

Figure 5. Schematic displaying the topological features of a stability diagram for reversible transitions. The two uniform arrangements are labelled by $U_1$ and $U_2$ and the intermediate line slip is labelled LS. Solid lines represent discontinuous transitions, continuous transition are marked as dashed lines and the arrows at the lines indicate the direction of the transition. For this type of transition $U_1$ transforms reversibly into $U_2$, showing only hysteresis above a threshold pressure $p_3$. (Taken from previous publication [29])

Below pressure $p_3$ the uniform structure $U_1$ can be transformed by a continuous (i.e. reversible) transition into its associated line slip. A further increase in D/d leads to a second continuous transition, whereby the line-slip is transformed into the uniform structure $U_2$. For a fixed pressure, the transitions occur at the same value of D/d, independent whether taking the forward trajectory (i.e. increasing D/d) or backward trajectory (i.e. decreasing D/d). This is indicated by the use of double-headed arrows in Figure 5.

For pressures p with $p_3<p<p_2$ the transition from the uniform structure $U_1$ to the line-slip remains continuous on both the forward and reverse trajectories. However, the transitions between the line-slip and the uniform structure $U_2$ are now discontinuous and no longer occur at the same value of D/d.

Increasing the pressure further (i.e. $p_2<p<p_1$) results in the transition from the uniform structure $U_1$ to the line-slip becoming discontinuous; furthermore the line-slip is no longer observed on the reverse trajectory (instead the system makes a discontinuous jump from $U_2$ to $U_1$). Finally, above $p_1$ the line slip also disappears from the forward trajectory and instead there are only discontinuous transitions between the two uniform structures $U_1$ and $U_2$.

To demonstrate the interpretation of Figure 5 we now give an example. Suppose we fix the pressure to be in the regime $p_3<p<p_2$ and start with the structure $U_1$ to the left of the diagram (i.e. low D/d). Keeping p constant and increasing D/d allows us to move towards the right of the diagram and eventually brings us into region C in the diagram. At the dashed boundary the change in enthalpy in going from $U_1$ to LS is continuous (as indicated by the dashed line). Increasing D/d (i.e. moving further towards the right) brings us into region B but the boundary can be ignored since the arrow points in the opposite direction. The final boundary (with an arrow pointing to the right) however is relevant and the solid line indicates a discontinuous change in the enthalpy, upon crossing the boundary the line-slip is transformed into the uniform packing $U_2$.

Reversing the trajectory we start with $U_2$ and decrease D/d (i.e. moving now towards the left), we first enter the region labelled B but the boundary can be ignored (since the arrow points in the opposite direction). Upon crossing from region B to C the solid boundary (with an arrow pointing to the left) is relevant and indicates a discontinuous transition from $U_2$ to the line-slip. Further decreasing D/d results in a continuous transition from the line-slip and back to the the uniform packing $U_1$.

The critical pressures $p_1,p_2,p_3$ depend on the choice of the uniform structures $U_1$ and $U_2$. and were determined from our simulations.

Figure 5, together with Table 1 containing the values for the critical pressure thus presents all the information required to characterise reversible transitions with $D/d<2.713$ (excluding structures with inner spheres).

5.2. The irreversible case – structural hysteresis

Irreversible transitions also lead from an initial to a final uniform structure (by a change of D/d), however, by definition the reverse path is not observed: it is unfavourable.

Within the range of cases discussed here we identify two types of irreversible transitions. The topological features of the first type, I, for which the initial uniform structure $U_1$ has the phyllotactic indices of the form $(l,l/2,l/2)$, is shown in Figure 6 (a). The second type, II, for which the initial uniform structure $U_2$ has phyllotactic indices of the form $(l,l,0)$ is shown in Figure 6 (b). Note that a common feature of both types of irreversible transitions is that the initial uniform structure is achiral.

Figure 6. Schematic stability diagrams for irreversible transitions. The uniform arrangements are labelled by $U_1$, $U_2$ and $U_3$ and the intervening line slip, in each case, is labelled LS. Solid lines represent discontinuous transitions, continuous transition are marked as dashed lines and the arrows at the lines indicate the direction of the transition. (a) Type I transition: Structure $U_1$ is of the type $(l,l/2,l/2)$. (b) Type II transition: Structure $U_2$ is of the type $(l,l,0)$.

5.2.1. Irreversible transition of type I

Starting with a uniform structure $U_1$ of the type $(l,l/2,l/2)$ we find that below $p_3$ the structure can be transformed continuously (and reversibly) into its associated line-slip, as shown in Figure 6(a) by the symbol ${\mathrm{U}}_1\leftrightarrow \mathrm{L}\mathrm{S}$. A further increase in D/d leads to the formation of a new contact, resulting in the uniform structure $U_2$. While this second transition is continuous (since the enthalpy varies smoothly throughout) it is not reversible, as indicated by the symbol $\mathrm{L}\mathrm{S}\to {\mathrm{U}}_2$

The transition from $U_1$ to $U_2$ is irreversible in the sense that beginning with $U_2$ a decrease in D/d does not lead back to $U_1$. For example, increasing D/d leads from $(4,2,2)$ to $(5,3,2)$ but a decrease in D/d from $(5,3,2)$ leads eventually to $(4,3,1)$ (and not back to $(4,2,2)$). This is evident also from the directed network shown in Figure 4).

For higher pressures in the range $p_3<p<p_1$ the transformation from the line-slip to $U_2$ becomes discontinuous (i.e. accompanied by a sudden change in the enthalpy). Finally, for pressures above $p_1$ the intervening line-slip between the two uniform structures disappears completely and this leads to a discontinuous change from $U_1$ to $U_2$.

5.2.2. Irreversible transition of type II

In the case of initial structures of the type $(l,l,0)$ the topological features of the transitions are similar, except that the transition is from the initial uniform structure $U_2$ to the final uniform structure $U_3$ and involves a decrease in D/d.

However, in the case of the $(l,l,0)$ structures a protracted difficulty was encountered in the regions indicated by ‘$\cdot \cdot \cdot$’ in Table 2. We expect the transition from $U_2$ to the line slip to be continuous below a critical pressure $p_3$, however we encountered pressures that were so low that we reached the limit of our numerical capabilities and were unable to obtain an estimates.

Table 2. Critical pressures for irreversible transitions.

Transition $U_1\Rightarrow U_2$	$p_1/{10}^2$	$p_3/{10}^2$	Transition $U_3\Leftarrow U_2$	$p_2/{10}^4$	$p_3/{10}^4$
			$(3,2,1)\Leftarrow (3,3,0)$	0.16	$\cdot \cdot \cdot$
$(4,2,2)\Rightarrow (5,3,2)$	2.89	1.62	$(4,3,1)\Leftarrow (4,4,0)$	0.62	$\cdot \cdot \cdot$
			$(5,4,1)\Leftarrow (5,5,0)$	0.16	$\cdot \cdot \cdot$

Note: The values for the pressures $p_1$, $p_2$ and $p_3$ determine the type of transitions encountered upon a change of the size ratio D/d in the stability diagram of Figure 6. This summarises the characteristics of all irreversible transitions for structures without inner spheres ($D/d<2.713$). The direction of the arrow indicates whether the transition happens at increasing D/d (marked by ⇒) or at decreasing D/d (marked by ⇐). The uncertainties of all pressures are again $\mathrm{\Delta }p=4\cdot {10}^{-5}$ where ‘$\cdot \cdot \cdot$’ indicate that the transition was not directly observed due to the numerical limits of our algorithm.

A further decrease in D/d eventually leads to a continuous (and reversible) transition from the line-slip to $U_1$. In the regime $p_3<p<p_2$ the transition from $U_2$ to the line-slip becomes discontinuous (as we have observed in simulations). Finally above $p_2$ the line slip disappears completely and the transition form $U_2$ to $U_3$ is by a single discontinuous jump in enthalpy.

6. Conclusion

When equal-volume hard spheres, bubbles or other similar entities are confined in narrow cylinders, they form a fascinating variety of helical and non-helical structures.

Extensive numerical simulations of such columnar packings showed that for monodisperse hard spheres the ratio of the cylinder diameter to the sphere diameter is the key parameter which determines the type of arrangement. For soft spheres both applied pressure and size ratio play a role and we have provided a complete taxonomy of the competing structures.

Here we have completed the survey of the hysteretic effects predicted for soft spheres, as was anticipated in our earlier work [29]. Transitions between packings, due to either changes in D/d or pressure P, are shown to be either reversible or irreversible. These relationships are summarised by a directed network and key parameters have been tabulated. Our results are found to be in broad agreement with other studies of transitions between columnar structures, either in the context of foams [6, 7] or colloids packed into tubes [14].

Our work is of relevance to the field of discrete or droplet microfluidics, in which ordered assemblies of bubbles or droplets are transported through channels [34–36]. Our directed network offers guidance how particular structures may respond to changes in either channel diameter or pressure.

A further relevance of our work may be to the recent observations of helical packings that form during the clogging of spherical particles in cylindrical channels [37].

Acknowledgments

We acknowledge inspiring discussions with M. Mungan.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

JW was recipient of an Irish Research Council Postgraduate Scholarship (Project ID No. GOIPG/2015/1998). SH was supported by the European Space Agency SciSpacE programme (project ‘REFOAM’, contract number 4000129502/20/NL/PG). AM acknowledges the support of the Supercomputing Wales project, which is part-funded by the European Regional Development Fund (ERDF) via Welsh Government.

Appendix

Simulation model

We consider the same model as in our previous work [29]. The interaction energy between spheres of diameter d is then given by, (A1) $E_{ij}^S=\{\begin{array}{cc}\frac{1}{2}(r_{ij}-d{)}^2& \text{if }{r}_{ij}\le d\\ 0& \text{if }{r}_{ij}>d\end{array}$(A1) where $r_{ij}=|{\mathbf{r}}_i-{\mathbf{r}}_j|$ is the distance between the centres of the spheres. Note the interaction energy is zero when there is no overlap between the spheres.

The interaction energy between the ith sphere and the cylindrical boundary is given by, $E_i^B=\{\begin{array}{cc}\frac{1}{2}k(r_{iB}-d/2{)}^2& \text{if }{r}_{iB}\le d/2\\ 0& \text{if }{r}_{iB}>d/2\end{array}$where $r_{iB}=|D/2-r_i|$. The spring constant k determines the softness of the spheres.

We conduct simulations using a unit cell of length L (volume $V=\pi (D/2{)}^{2}L$), containing N spheres. On the ends of the unit cell we impose twisted periodic boundary conditions. The periodic boundaries are implemented by placing image spheres above and below the unit cell, where each sphere of the unit cell is moved in the z-direction by L (and $-L$) and twisted by an angle α (and $-\alpha$, respectively) in the xy-plane [27].

Our simulations were run with 12 spheres in the unit cell, commensurable with the minimum number of spheres required for the previously identified hard-sphere packings in the investigated range of D/d (N = 1, 2, 3, 4, see Table 1 in [27]). For simulations of the $(5,5,0)$ structure or its adjacent line slips, which are the only structures with at least five spheres in the unit cell, N = 60 spheres was used.

Stable structures can be generated by minimising the enthalpy H = E + PV. Here $E=E^S+E^B$ is the internal energy due to overlaps with other spheres ($E^S=\sum _{i\le j}^N{E}_{ij}^S$) and with the cylindrical boundary ($E^B=\sum _i^N{E}_i^B$)) and P is the pressure. During the minimisation, the free parameters are the sphere centres $\left\{{\overrightarrow{r}}_i\right\}$, the twist angle α, and the length of the unit cell L, while the pressure P is kept constant.

Enthalpy and pressure have to be rescaled accordingly to obtain non-dimensional quantities. We use the dimensionless enthalpy $h=H/\left(k{d}^2\right)$ and dimensionless pressure $p=P/\left(k/d\right)$, where k is the spring constant and d is the sphere diameter.

From our previous results we are able to identify the minimal enthalpy structure for a given value of p and D/d. For this stable structure we ask how does it change when the pressure p and diameter ratio D/d are continuously varied (where in the latter case we only explore the change in D/d by varying D and not by varying d), while remaining in the local minimum of enthalpy? Thus, we will only be using the BFGS method [38] here, which uses a conjugate gradient algorithm to search for the nearest local minimum.

With this procedure, structures of high symmetry, such as the $(4,2,2)$ or $(3,3,0)$, can get stuck on a saddle point, where the enthalpy is not minimal, but its gradient is zero. This can be avoided by applying a small random perturbation to the structure, which displaces it from the saddle point, followed by a local minimisation. Starting with a given uniform structure, a change in D/d may lead to a mechanically unstable structure which is identified as a structure that starts out as the expected line-slip arrangement but then collapses into a metastable structure.