Renaissance of Bernal's random close packing and hypercritical line in the theory of liquids

We review the scientific history of random close packing (RCP) of equal spheres, advocated by J D Bernal as a more plausible alternative to the non-ideal gas or imperfect crystal as a structural model of simple liquids. After decades of neglect, computer experiments are revealing a central role for RCP in the theory of liquids. These demonstrate that the RCP amorphous state of hard spheres can be well defined, is reproducible, and has the thermodynamic status of a metastable ground state. Further evidence from simulations of square-well model liquids indicates an extended role of RCP as an amorphous ground state that terminates a supercooled liquid coexistence line, suggesting likewise for real liquids. A phase diagram involving percolation boundaries has been proposed in which there is no merging of liquid and gas phases, and no critical singularity as assumed by van der Waals. Rather, the liquid phase continuously spans all temperatures, but above a critical dividing line on the Gibbs density surface, it is bounded by a percolation transition and separated from the gas phase by a colloidal supercritical mesophase. The colloidal-like inversion in the mesophase as it changes from gas-in-liquid to liquid-in-gas can be identified with the hypercritical line of Bernal. We therefore argue that the statistical theory of simple liquids should start from the RCP reference state rather than the ideal gas. Future experimental priorities are to (i) find evidence for an amorphous ground state in real supercooled liquids, (ii) explore the microscopic structures of the supercritical mesophase, and (iii) determine how these change from gas to liquid, especially across Bernal's hypercritical line. The theoretical priority is a statistical geometrical theory of RCP. Only then might we explain the coincident values of the RCP packing fraction with Buffon's constant, and the RCP residual entropy with Boltzmann's ideal gas constant.


Introduction
A hitherto unattained objective of liquid-state theory is the development of a formalism that can predict the thermodynamic state functions and liquid phase bounds from a model Hamiltonian defined by intermolecular potentials. In essence, the problem is one of predicting the relative positions of molecules of a given liquid at a given (p,T ) state point in a way that enables thermodynamic properties that depend upon this structure to be described and computed analytically.
In an ideal crystal, the positions of molecules in a unit cell can be specified precisely, and the structure of the extended crystal reproduced simply by repeating that unit cell in three dimensions. The structure of a real crystal can be approached through perturbations from this ideal, for example thermal motions and the presence of structural defects such as vacancies and dislocations. The underlying periodicity simplifies the calculation of the crystal properties from known molecular parameters.
An ideal gas is in the opposite extreme; the instantaneous positions of the constituent molecules are totally random. Just as we approach a real crystal by perturbing the regularity of the ideal lattice, we approach the non-ideal gas by describing the spatial correlations, i.e. by increasing the order. When the molecules are given finite size and intermolecular interactions, there is a departure from the perfect randomness of the structure as pairs of molecules interact. A well-known approximation is the van der Waals equation [1]. As the density of the gas is increased further, we need to consider the effects of higherorder (triplet, quadruplet, . . . etc) correlations, but so long as the density remains reasonably low, a non-ideal gas cluster expansion predicts to good approximation the structure and properties of a simple gas, e.g. argon, composed of weakly interacting atoms [2].
It is not surprising, therefore, that the early theoretical attempts to develop a theory of liquids approached the problem from both of the above extremes: by considering the liquid to be either a highly disordered crystal [3] or a high-density gas [4]. Noting that a simple liquid close to its melting point has a density of the order 90% of its corresponding crystal, it seems reasonable that a theory based on a disordered crystal can reproduce liquid properties such as the density and energy. Likewise, treating a liquid as a dense gas can give reasonable predictions of its entropy. But the gas approach fails where the crystal approach succeeds-and vice versa. Crystal-based models, because of the regularity of the underlying crystal lattice, have far too little entropy. And cluster expansions [2], which describe non-ideal gases, fail to converge as they approach the density of the liquid and cannot therefore successfully describe dense liquids [5].
J D Bernal, one of the great polymaths of the 20th century, who was responsible for much of the early development of crystallography, had a core interest in the biological role of water. He argued that if he was to understand how water was relevant in biological processes he needed to understand the structure of water. Together with R A Fowler in 1933, Bernal published an extensive paper on the structure of water and aqueous solutions [6]. Realizing, however, that water was a difficult problem, he subsequently commented [7]: 'It is not worth tackling complicated liquids until we understand simple ones'.
Bernal therefore shifted his focus to simple liquids exemplified by liquid argon. The first of two attempts [7,8] he made at the problem in the 1930s [8] was based on an analogy made by Zernike and Prins [9], who pointed out that there was similarity between the x-ray scattering pattern of a liquid (smooth but with broad peaks) and that of a crystalline solid (sharp peaks). Bernal interpreted this analogy [8] as 'the peaks of the first, second and third coordination spheres corresponding to the reflexions from planes of the same intermolecular distances as in a crystal'. Bernal found this to be a delusive approach because it overestimated the degree of long-range order. In his second approach, published in 1937 [7,10], he tried to describe a liquid in terms of the average co-ordination of the atoms. In retrospect [8], once he realized that the values required could not be calculated theoretically, he recognized this as a 'highly artificial' unquantifiable description.
After a gap of about two decades, in the 1950s Bernal returned again to the problem. After hearing a lecture by F C Frank on atomic packings in complex metal alloys [11], he began to develop his concept of liquids as 'homogeneous, coherent and essentially irregular assemblages of molecules containing no crystalline regions nor, in their low temperature form, holes large enough to admit another molecule' [7,8]. As will be described below, this concept was ultimately physically realized as the random close packing (RCP) of equal spheres.
This last statement of Bernal is especially pertinent to our current understanding of liquids, which we describe in later sections. The probability of finding 'holes' in a liquid is in fact a property of central theoretical importance. It is its chemical potential which determines the thermodynamic state of existence [12][13][14]. We will argue that the liquid phase ranges from a metastable RCP state at 0 K, or infinite pressure for hard spheres, with zero holes, to a percolation transition at a lower density whereupon holes coalesce to become one large hole, or 'nano-bubble of gas' that marks the end of the liquid phase. Moreover, we will see towards the end of this review that as a liquid approaches its triple point, even at equilibrium there is evidence of fluctuating heterogeneities or 'nucleites' that were postulated over half a century ago. Not even Bernal could foresee the heterophase fluctuations, resembling incipient gas and molecular microaggregates, that now appear essential aspects of liquids, especially in the vicinity of phase bounds.

Early models of RCP
The earliest attempt to realize a random packing model was that of Westman and Hugill [15] in relation to ceramics, followed by Prins [16] in relation to liquids. Prins dropped seeds onto a glass plate, photographed them and by counting the number of seeds in concentric rings from central seeds constructed a radial distribution function (RDF). This is the quantity that is primarily used to describe the structure of a liquid, and is defined as the probability of an atom being at a given distance from a central atom. Although limited to two dimensions, the RDF obtained resembled the experimental result for liquid mercury. Menke performed similar 2D experiments using steel spheres [17], and obtained an RDF similar to that that he had measured for liquid mercury. Morrell and Hildebrand were perhaps the first to tackle the problem in three dimensions [18], when they used coloured gelatin spheres to represent atoms, cleverly counteracting the influence of gravity by using a gelatin solution as the surrounding medium.
None of these early experimentally produced packings achieved densities close to today's generally accepted packing density of RCP close to 0.64. Later work by O K Rice [19] obtained the density of a RCP state through measuring the volume of water required to fill the interstices of a packed bed of glass spheres. After making corrections for the perturbing effects of the container, he concluded that his random packing had a volume some 15(±3)% greater than that of crystalline cubic close packing. This value matched the expected density change on melting of an inert gas crystal.
By the 1960s computers had begun to make simulation of random packings possible, although limited by the available computational power to relatively small systems. In 1955 Alder et al [20] performed Monte Carlo (MC) calculations on periodic systems of 80-100 hard spheres. Although he had no direct experience of using digital computing, Bernal was, however, able to see the possibilities of computer simulations. He founded a computing laboratory and used its capabilities to simulate the procedure by which he had built a realization of his concept of a random packing in the laboratory [21]. The co-ordinates of points were chosen randomly; those that were greater than a specified distance from the points already chosen were retained [7]. Bernal noted that as the method proceeded, it became very inefficient-hundreds of points had to be chosen in order for one to be accepted. The result of this process is what we now call 'maximum random parking' [22]. This non-equilibrium state of hard-spheres becomes, de facto, a reference point for non-crystalline close-packed states [23], a range of states of which Bernal's RCP is now seen to be the limiting equilibrium state, albeit metastable with respect to crystal phases.

Laboratory and computer realizations of RCP
Bernal's initial attempt to construct RCP shook together several thousand 1/4 ball bearings. This assembly was further densified by winding thick rubber bands around it and then kneading the mass [24]. In addition to obtaining the density of maximum RCP, he also characterized a lower density state, obtained without shaking and compressing, which he termed 'random loose packing'. From these assemblies Bernal obtained information on local coordination geometries, noting in particular that the distribution in number of contacts was distinct from that of either of the regular close-packed crystalline arrangements (face-centred and hexagonal cubic close packings), where every ball has twelve contacts. Variation of contact numbers was found to be a significant feature of the irregular random-packed noncrystalline arrangement.
Bernal made comparisons between the properties of his conceptually simple model and those of liquid argon [24]. If the RCP structure really was a plausible model of the structure of the liquid at its melting point, the volume increase on melting should be consistent with experimental data on argon. It was. From these early results, he also estimated the latent heat of melting to be between 1/4 and 1/6 of the energy of evaporation; the observed ratio is actually 2/11. Considering random loose packing as a model of an inert gas liquid at a higher temperature, his calculation of the difference in configurational energy between the liquid at its melting and critical points was within about 12% of the experimental value.
Experimental random packings depend on container size and shape; structures are more ordered in the vicinity of the bounding surfaces. G D Scott, who also became involved in random packing studies in the 1960s, solved this problem [25] by taking density measurements of packings in vessels of different sizes; by extrapolating to infinite size he concluded that RCP had an upper density limit of 0.637, while 'random loose packing' had a lower stability limit of 0.601. Later work on larger packings from both laboratories confirmed this upper density limit as 0.6366 ± 0.0004 [26,27] or ±0.0008 [28]). It was observed that this value is tantalizingly close to Buffon's constant 2/π (0.6366192) [29], a constant that occurs naturally in solutions to problems in statistical geometry.
Bernal and Scott [30] collaborated in calculating RDFs of their models using both sets of data. The earliest comparison they made between RCP and experimental results from neutron scattering for liquid argon is reproduced in figure 1. The agreement, both between the two independently constructed laboratory models and with the experimental data, was quite striking at the time.
This comparison demonstrated that Bernal's hypothesis of 'homogeneous, coherent and essentially irregular assemblages' was indeed consistent with experiment, and so answered his foremost question: 'what is the essential structure of a simple liquid?'.
When numerical coordinate data were measured in subsequent work, not only could more quantitative comparisons be made with experiment, but also the influence of the softness of the potential functions between 'real' atoms could be explored. For example, by assigning appropriate Lennard-Jones potentials to the innermost spheres of Scott's model, Bernal and Finney obtained good estimates of not only the volume changes on melting but also the heats of fusion of neon, argon, krypton and xenon at pressures up to 6000 atmospheres [31,32].
In order to obtain better statistics for the RCP structure, Bernal wanted a larger sphere assembly with more-precisely measured coordinates. This was realized in the 'final' large (7934 sphere) random close-packed model that was built in Bernal's laboratory [26,27,33]. The most important conclusions drawn from that model were [26,27] (a) its density (0.6366 ± 0.0004 mentioned above); (b) its RDF and (c) its statistical geometrical structure. The RDF (figure 2) clearly shows a split second peak, a feature that earlier 'low resolution' models were unable to show. This splitting is not observed in real liquids as it is smoothed out by the softness of their intermolecular potentials-though it is observed in many metallic alloy glasses.
This large model provided the data for developing a statistical geometrical characterization of random packing in terms of the topologies of the equivalent Voronoi polyhedron assembly 3 . Properties that were characterized included the average and distribution of number of faces per polyhedron and the average and distribution of edges per polyhedron face. At a more detailed level, polyhedron types F n were defined by the number F of faces of n sides, or [F 3 ,F 4 ,F 5 ,F 6 ,F 7 ,]. Using this descriptor, polyhedron type distributions for the Scott model were found to be similar to those of Bernal's large model [27], implying that the two models had the same structure. The Voronoi characterization has also been useful in demonstrating differences between the structures of differently constructed computer models [34]. Estimates of the configurational entropy of RCP using polyhedral type as a state definition were only 3% lower than the experimental value for liquid argon [35].
This work on laboratory-built models was performed largely in the early to mid-1960s when computer resources were insufficient to construct dense packing models by simulation.
The computational resource needed to accommodate the large number of local, co-operative structural rearrangements that would simulate the physical shaking or compression that is relatively easy for physical packings was not then possible. It was, however, feasible to construct models by sequential deposition of spheres. One of the earliest to explore such techniques was Bennett [36] who obtained RDFs similar to those of Bernal's large laboratory model. Moreover, although computational limitations restricted them to constructing finite clusters only, Mason [37] and Finney [38] independently implemented a hard-sphere gas compression procedure that resulted in models with RDFs that were consistent with that of the large laboratory model.
With increased computer power, it became feasible to simulate the construction of dense packings using MC and molecular dynamics (MD) methods of computational statistical mechanics. The earliest computer simulation of RCP using periodic boundary conditions, a device that is not available to the laboratory modeller, generated a spacefilling packing of density 0.6372 ± 0.008, within estimated uncertainty indistinguishable from the 0.6366 ± 0.0004 of the large laboratory model [39,40]. The RDF was also consistent with that of the large laboratory model, showing no peak at √ 2 diameters (a characteristic distance in a close-packed crystal) and, in particular, the splitting of the second peak into two peaks close to √ 3 and 2.

Liquid state theory and RCP
The construction and characterization of the large laboratory model took place just as the computer simulation of simple liquids was becoming well-established technology. The detailed microstructure of any model Hamiltonian of a simple liquid could then in principle be obtained and analysed. Using the Voronoi characterization tools, the local structure of the RCP state could now be compared with that of a computer simulation of a simple liquid. Additionally, the structural consequences of the impenetrability of the hard-sphere potential of the RCP model could be explored. Consequently, simulations were undertaken of liquid argon just above and below its melting temperature, together with other computations using potentials with increasingly hard repulsive cores [41]. The results of these calculations indicated an RCP origin of the liquid phase structures. In figure 3 we show a present-day comparison of the RDF of a high-resolution computer-generated RCP state [42], and compare both with the RDF of the Lennard-Jones model of liquid argon at its triple point. The first observation is the similarity of the computer-generated RDF with the RDF of Bernal's model shown in figure 2. The RCP state is characterized by a complete absence of a peak at r/r o = √ 2 (a characteristic close-packed crystal distance) and the existence of peaks close to √ 3 and 2. In the equilibrium triple-point liquid there is just a single second peak, but on quenching the liquid to a glass at almost 0 K (T * = k B T /ε = 0.05) and then annealing, the second peak splits as seen in the RCP state. The two states belong to the same 'liquid phase'. This modern comparison of Bernal's RCP with simple liquids vindicates his original contention. We will show in the following sections that the RCP of hard spheres is indeed the appropriate reference state for a modern theory of simple liquids.
In his treatise on the physics of disorder [43], distinguished theoretical physicist John Ziman states 'This simple idea is now seen to be the key to any qualitative or quantitative understanding of the physics of liquids'. A complementary observation was made in 1970 by Rowlinson in his introductory Faraday Discussion lecture [44]: 'It has therefore been hard to admit that the form, or even the existence, of the attractive forces has little direct effect on the structure of a liquid... The recent realization of this truth has followed the extensive studies . . . of the properties of assemblies of hard spheres without attractive forces'. Subsequent successful perturbation theories [4,45] further underline the veracity of Rowlinson's comment.
All of the above historical observations raised the hitherto neglected question: what is the thermodynamic 4 status, if any, of RCP? Furthermore, what is the significance of 'random loose packing' and of 'maximum jammed states?' Recent computational research indicates that the state of RCP has in fact a well-defined thermodynamic status. This has profound implications in the general theory of liquids and glasses. It is to this issue we now turn; we will argue that all other 'jammed' 4 Throughout this review we use the adjective 'thermodynamic' to describe the state of a molecular system, or model Hamiltonian, that is thermally equilibrated and can be represented as a state-functional point on a Gibbs T , p surface, e.g. of density or entropy. This implies that there can be no gradients of T , p or ρ, and that the system is in a state of equilibrium or metastable (local) equilibrium whereupon its Gibbs energy G T ,p is a minimum, and for any process of change to any other local state of the system dG T ,p > 0. states, excepting crystalline structures, by contrast, have no thermodynamic status.

Excluded and accessible volume
The simplest imaginable model fluid of rigid elastic spheres has played a central role in the development of liquid state theory since the pioneering work of van der Waals [1]. He obtained an equation-of-state resembling that of real fluids by simply adding an attractive potential to a hard-sphere fluid equation which he obtained by simply assuming an ideal gas volume could be replaced by the excluded volume on collision of two spheres: where Z is pV /(N k B T ), ρ is the number density Nσ 3 /V and b 2 = 2πσ 3 /3. b 2 is the excluded volume (V e /N ) and is the second virial coefficient of the hard-sphere fluid. The complement of V e is the available volume (V -Nb 2 ). When this equation is expanded in powers of density, all the virial coefficients take the value b 2 (equation (1)). However, a formally exact virial expansion of gases obtained from the cluster theory of Mayer [2]: shows that only the first two terms of equation (1) are correct. In fact, equation (1) is inadequate except at low gas densities. Virial coefficients beyond n = 2 can be computed by integrating the many-body excluded volumes of clusters of n spheres over all space. Implicit in current mainstream theories of simple liquids [4], moreover, is the assumption that equation (2), given all the virial coefficients, represents the hard-sphere fluid thermodynamic pressure in the density range of real liquids up to freezing. An alternative expression for the accessible (or excluded) volume of equilibrium hard-core fluids that relates directly to the chemical potential was derived by Hoover and Poirier [12]. In the thermodynamic limit of a very large system, V a can be defined as the volume accessible for the relocation of any one sphere. This is exactly related to the excess chemical potential (µ) by the equation of Hoover and Poirier [12], which also relates to the Mayer cluster expansion of the chemical potential in powers of density: The summation in equation (3) is only known to be exact for the hard-sphere fluid pressure in the low density limit [2]. Modern treatises of the theory of simple liquids, however, invariably begin with the assumption that equation (3) represents the hardsphere reference fluid up to high densities that correspond to the densities of real simple liquids [4]; this assumption is now seen to be invalid.

Percolation and percolation transitions
If we continuously reduce the volume occupied by a number of entities-for example molecules of a gas-there comes a point at which sufficient contacts between the entities have been made for a continuous path between the entities to span the system. This system is then said to percolate, and the state point at which this path appears is a percolation transition. The physics of percolation transitions is much simpler and easier to access for discs than for spheres For 25 years there have been postulates of thermodynamic discontinuities at percolation transitions [46][47][48] but the formal investigation of these properties for hard-sphere fluids is relatively recent [49,50]. The first evidence for identifying a percolation transition with a thermodynamic phase transition was obtained from an MD study of hard parallel cubes [47]. The available-volume percolation transition-the point at which a continuous path exists through the available volumewas found to coincide with a weak discontinuity in transport coefficients as functions of state, and a deviation of the virial equation-of-state from the thermodynamic pressure. Using scaling arguments, Kratky [48] was able to make rough estimates of hard-sphere percolation transitions. These he discussed in a paper with the searching title 'Is the percolation transition of hard spheres a thermodynamic phase transition?'.
Values for the densities of the available-volume (V a ) and of the excluded-volume (V e ) percolation transitions have been reported for the hard-sphere fluid [48,50]. V a is the volume available for insertion of a sphere whereas the excluded volume V e = V − V a , i.e. the volume of an excluded zone of radius 1σ from each sphere centre. This is simply illustrated for hard discs in figure 4. In two-dimensional fluids, the percolation transition density of V e coincides with the percolation transition density of V a . Thus the excluded volume and available volume percolation densities y pe and y pa are equal, i.e. for D = 2, y pe = y pa = 0.31 [48]. This is not the case in three dimensions where the two percolation transitions in the hard-sphere fluid occur at different packing densities y pe = 0.046 and y pa = 0.281 [50].
The breakthrough towards answering Kratky's question [48] for the hard-sphere fluid came with the application of high-performance computing to determine virial coefficients b 8 to b 10 [51]. This has enabled the hard-sphere equation, for both discs and spheres [52], to be written in closed form with very high precision [50]:

of equation (3) by Clisby and McCoy in 2006
where Z = pV /NkT and ρ * = ρ/ρ 0 . The constant is empirical to within four figures. Equation (4) with m = 8 reproduces the Padé approximant to the full virial series [51] over the full density range. However, it reproduces the thermodynamic pressure only up to the density of the available volume percolation transition (y pa = 0.281) with 6-figure accuracy. This observation has profound implications for the development of a general theory of liquids.
These highly accurate virial equations have stimulated the calculation of equally accurate, essentially 6-figure : the left hand side shows a high-density liquid-like region above the percolation threshold density; the right hand region is gas-like below the percolation threshold density. Reproduced from [91] with permission of The Royal Society of Chemistry. precision, thermodynamic data by MD computation. Using the high-performance hard-sphere program DYNAMO [53] a bifurcation of the pressure given by the Mayer virial equation and the thermodynamic pressure has been observed. When the virial pressure is compared with the thermodynamic pressure, a deviation is discernable near to the available volume percolation transition (PA) [54]. In figure 5 we show the virial equation-of-state and the Bannerman DYNAMO data up to freezing density. The virial expansion is essentially exact, albeit numerical, with 6-figure accuracy at low density. Convergence up to the reduced density ρσ 3 = 0.025 has been rigorously proven [55]. Thus, when the thermodynamic MD pressure deviates from the virial equation, there must be a thermodynamic phase transition because the virial series is mathematically continuous up to a pole at crystal packing (equation (4)). Figure 5 shows that there is a bifurcation at a density no greater than y pa .
Given that there is a thermodynamic discontinuity at a density well below fluid freezing, how does this relate to liquid-state theory reference states? Historically [1][2][3][4], from van der Waals to modern mean-field theories, and also including integral equation approaches such as Percus Yevick and CHNC [4], there has always been a central hypothesis: these theories all assume that the virial expansion is continuous to liquid densities. Figure 5 shows that this assumption is not valid; the virial equation of gases fails to converge on to the thermodynamic equation-of-state at the available volume percolation transition density (y pa ). Beyond this point an alternative approach to the equation-of-state of the dense hardsphere fluid, and hence to the liquid state also, is required.
It is a property of thermodynamic state functions on a Gibbs density surface, such as figure 5, that the density will be a continuous state function in all its derivatives. Moreover, the freezing transition is a first-order thermodynamic phase transition and as such, each coexisting phase has a theoretical extrapolated existence on either side of the phase transition. These metastable thermodynamic states exist, but may not be accessible experimentally. Nevertheless, we can deduce that a continuous state of the hard-sphere fluid, from equilibrium liquid-like structure for densities above the percolation transition, to metastable states beyond the equilibrium freezing transition, to an RCP ground state at T * → 0, belongs to the same thermodynamic phase-and hence same continuous equation-of-state-as the liquid-like region of the hard-sphere fluid as shown in figure 5.
Consequently, there exists an idealized reference state for the development of the theory of liquids and the thermodynamic equation-of-state of liquids. This RCP state is the liquid-state analogue of the ideal gas and the perfect crystal, which themselves are idealized states for theories of the gaseous and crystalline phases respectively. There is thermodynamic continuity of the gas phase up to the density of the percolation transition that bounds its existence. So likewise for the crystal phase with continuity from the perfect lattice up to a melting line. A consequence of the thermodynamic continuity, from equilibrium dense fluid to metastable RCP, is that the starting point for theories of simple liquids is not the ideal gas, as has hitherto been generally assumed [1,2]: it is an 'ideal glass'.
Next, we investigate the evidence that this limiting metastable amorphous ground state of the hard-sphere fluid is indeed the same RCP structure that was originally investigated by Bernal.

Processes leading to RCP
We have seen from the introductory sections and figures 1-3 that RCP has long since been well characterized, but is it well defined? If so, what is its thermodynamic status? Any state of matter not in thermodynamic equilibrium can only be as well defined as the process that produces it. To ascertain the answer to this question, we examine in the following section specific processes, both irreversible and reversible, that lead to RCP structures.

Irreversible
There is an extensive literature describing various geometrical constructions of random close-packed structures by both analytical and numerical mathematical methods [56,57]. None of these mathematical procedures, however, contain evidence of thermodynamic status. Such evidence would be proof of a local density or entropy maximum. Nonequilibrium, or metastable, solid states can be defined as thermodynamic states only by fixing two macroscopic state variables, and one or more process constants. For a phase in equilibrium, temperature and pressure are usually used to define the state. The thermodynamic properties such as enthalpy or entropy are state functions with a unique value. Most composites such as alloys, polymers, ceramics, etc., appear to be thermodynamically stable, but they will however have thermodynamic properties that depend on the process of production, and in many cases on the timescale of the process. Non-equilibrium solid states can only be well defined if the compaction or cooling process is itself well defined. Moreover, the initial states from which a compaction or cooling process begins must also be thermodynamic, i.e. well defined.
The earliest compaction computer algorithm [39,40] simply modified the MD equation-of-motion for hard spheres by fixing the rate at which the volume (V ) containing N spheres decreased by specifying an equivalent rate increase in (dσ 3 /dt) T where sigma is a dynamical sphere diameter at constant total volume and temperature (T ). Jodrey and Tory [58] devised a similar well-defined algorithm by a single rate constant, albeit athermal, starting from a dilute hard-sphere gas. Unfortunately, the scope of compression rates accessible at the time restricted their observation range of maximum packing densities from 0.64 to 0.65. An inspection of figure 6 of the paper by Jodrey and Tory suggests a slightly faster rate would have resulted in a 'thermodynamic' RCP state accessible as found previously [39,40] with this type of algorithm.
In many compaction processes, the rate of change of volume (dV /dt) T is kept constant. Such a process can be said to be zeroth order in 'free volume', whereupon the relaxation processes allow equilibration at low density, but rapidly freeze at high density. Whilst prohibiting crystal nucleation, such a process also misses the 'window' of re-equilibration to access RCP. There is a growing literature [58][59][60] of similar such processes. Some of these 'simulation' phenomena may have no experimental counterpart, and exist only within a computer. These algorithms result in closepacked structures that have been called 'maximum jammed states'. Here, the compaction rate is too fast to avoid nucleation; it is also too fast for re-equilibration onto the RCP state, and homogeneous amorphous 'maximum jammed states' for packing densities y > y RCP are obtained. Thus, this area of science may only exist in the computer algorithms without any obvious relevance to any real molecular level counterpart, although the jammed states can be identified with experimental properties of granular media [61,62]. If the compression rate is not high enough to avoid incipient nucleation, a glass transition may occur whereupon a range of heterogeneous glassy states is obtained. These processes resemble glass transition phenomena of real liquids [40]. We can define an isothermal compaction process by a rate law that is first-order in free volume [42]: where V is the sphere-containing cube volume that is subject to periodic boundary conditions. V ∞ is the close-packed state atkT/pσ 3 = 0. Though V ∞ is a priori unknown, at any point in time of the compaction process it is accurately predictable from the temporal pressure using a self-consistent free-volume equation of state which is very accurate at high density, and exact near the close-packed limit: where Z = pV /N k B T . The angular brackets denote a shorttime ( t) average. At constant total volume, the spheres are expanded at a constant rate dσ/dt, which is obtained from V (t) and V ∞ by This process of compaction is a simple modification of the collision-predictor algorithm for classical motion [40,42]; the quadratic equation for a time t to next collision becomes Equations (5)-(8) describe a well-defined compaction process; the final compacted state is thus defined by the firstorder compaction rate constant k 1 that defines (dσ/dt) in equation (7). The final state may require ensemble averaging of the starting state for small systems. For sufficiently large systems a single compaction of an equilibrium starting configuration suffices. Typical results [63] for the dependence of the limiting maximum density on compaction rates (log k 1 ) ranging from −6 to +2 are shown in figure 6. Limiting packing densities were determined for each compression simulation by extrapolating to Z = ∞ using equation (6). At very low k 1 the system can crystallize to the perfect face-centred cubic crystal for small periodic systems, where for the close-packed crystal density y 0 = 2 1/2 πNσ 3 /6 (0.7405). As k 1 increases we obtain defective crystalline states with reduced densities. There follows a range of steeply decreasing density, which probably denote heterogeneous-probably glassy-states with small crystal nuclei. Eventually an RCP plateau is observed over a limited range of k 1 . When correction was made for finite size, a value of 0.636 ± 0.002 was reported for the packing fraction of the RCP plateau, consistent with the accepted value for RCP.
Because the 1st order isotropic compaction process is completely homogeneous all configurations are 'force balanced', but in the non-equilibrium states there is 'mechanical instability' with respect to decrease in the configurational pressure. All the state points in figure 6 are, in a sense, 'jammed', simply by definition of infinite pV/NkT or zero T , i.e. approaching 0 K they become thermally static. For the low-density structures the rate of compaction is so fast that there is no time for reorganization. These structures are mechanically unstable with respect to relaxation when the compaction process is suddenly stopped. For higher values of k 1 , (order) 10-100, 'jammed' loose-packed amorphous states can be obtained. Such structures resulting from essentially athermal processes have no thermodynamic description. Moreover, they may not meet the same criteria for 'jammed' states, as defined by different compaction processes [58][59][60][61][62]. Beyond around k 1 = 100 the algorithm equations (5)-(8) becomes unstable near close packing when many collision times in equation (8) become less than the smallest precision of the compiler.

Reversible
It has also been demonstrated that an RCP state of the hardsphere fluid which is essentially the same as that obtained by irreversible compaction can be produced reversibly. This thermodynamic route to avoid crystallization makes use of a single-occupancy (SO) constraint that was originally devised by Hoover and Ree [64] to achieve just the opposite: the maintenance of a crystalline order in their computation of the freezing transition of hard spheres. The centre of each sphere is confined within the fcc primary Voronoi polyhedra.
A similar constraint (simple-cubic single-occupancy (SCSO)) has also been found to lead from a low-density constrained gas to an RCP state [63]. The hard-sphere simple-cubic (SC) crystal structure is thermodynamically unstable with respect to other crystal structures and also mechanically unstable with respect to slip into local more random arrangements. The SCSO fluid compacts reversibly from the SO-ideal gas at y = 0 to an RCP structure at y = 0.637, as shown in figure 7 All the data points in figure 7 are in thermodynamic equilibrium, though there is the possibility of more ordered structures that are not accessible in the MD computations. The SCSO equation of state is everywhere very close to the hard-sphere fluid, being essentially the same in the limits of both low density, where it obeys the ideal gas equation, and high density, where it follows the RCP free volume equation. Beginning with a low density, a repeated sequence of compression and re-equilibration leads to stable pressures for millions of collisions per sphere.
An important property of SO systems is the particlewall to particle-particle collision frequency ratio. At RCP it approaches a constant value of about 0.03. The total number of particle-particle collisions per unit time diverges with the pressure, but the ratio of particle-particle to particlewall collisions stays constant. It appears that the RCP state is being stabilized in the vicinity of close packing. On removal of the SO constraint in the immediate vicinity of close packing, there is an initial relaxation to a slightly lower pressure and consequently a slightly higher limiting RCP density is obtained. For a system of 64 000 spheres, a value of y RCP = 0.6365 ± 0.001 was reported [42,63], again the same as previous experimental RCP packing fractions from various sources.

RCP density and Buffon's constant
Some papers on RCP have asked, rather than answered, a number of questions about the definition and characterization of RCP [65][66][67]. To address the first of these questions: 'What is RCP [65]?', we first rephrase it as 'does RCP have a thermodynamic status?' From the previous section, we can say yes.
In response to 'Is RCP well-defined [66]?', we can also say yes. The RCP ground state of the hard-sphere fluid appears to be unique, but proof of such thermodynamic metastable status remains elusive, as also does any statistical theory of its structure. However, we do now, in principle, have another thermodynamic definition given by a limiting constant in a free-volume or cell-cluster expansion equation-of-state of the equilibrium hard-sphere fluid. We can in principle obtain it to any requisite accuracy in the amorphous high-density region and it connects apparently continuously to a RCP amorphous ground state.
We can now also address the question whether RCP produced in the above processes is a maximum or minimum in density and entropy. The definition of thermodynamic equilibrium, or metastable equilibrium, implies that all other possible states in the vicinity of RCP must increase the Helmoltz energy: RCP therefore represents a local minimum in free energy, which from equation (9) requires a local entropy and/or density maxima at constant temperature. This thermodynamic analysis can answer the third question: 'Why is RCP reproducible?' [67]. All thermal systems not in a state of equilibrium will tend towards the nearest local equilibrium. If there is a mechanism and provided the timescale for nucleation is avoided, any systems of spheres given some kinetic energy and time, will 'gravitate' towards the maximum in density and entropy. In any well-defined process however, with equipartition of energy prevailing, RCP may only be obtained if crystal nucleation is completely avoided.
A further unsolved query is the closeness of the value of RCP density of 0.6366 ± 0.001 to Buffon's constant 2/π (=0.6366192), which often crops up in analytical solutions of simpler solvable problems in the general field of statistical geometry [29].
This suggestion of a possible relationship to Buffon's constant arose from the earliest sphere packing experiments of Bernal and Mason [24] and Scott [25]. Experimenting with 20 000 steel balls contained within copper cylinders, Scott obtained a maximum density value of 0.634. This was later refined by Scott and Kilgour [28] to 0.6366, a result essentially identical to that obtained earlier from the large hard sphere model by Finney [26,27]. Another three values from more recent analyses have been added to the list. Using more recent MD pressure data, Le Fevre's interpolation to 1/Z = 0 gives y RCP = 0.6375 [63,68]. The first-order compaction process estimated for N → ∞ as discussed in section 6.1 gives y RCP = 0.6360 [63] ( figure 6). Finally, a value of 0.6365 is produced by relaxation of RCP-SCSO to an unconstrained RCP state ( figure 7).
There are many other 'RCP' computer-generated values in the literature. Using a wide range of algorithms, various jammed RCP states with packing densities from y = 0.60 to 0.66 have been obtained. Such states, however, are ill defined in a thermodynamic sense. In contrast, as we have seen, the thermodynamic RCP amorphous ground state is well characterized. It can be also well defined. However, there is at present no analytic theory of RCP. A consensus packing density from all sources suggests a value of 0.6366 ± 0.0005. Within the numerical uncertainties, this is essentially Buffon's constant. In the absence of any analytical theory of RCP, we wonder if it might be possible to construct RCP mathematically in such a way as to obtain y RCP = 2/π.

RCP entropy and Boltzmann's constant
A difference in entropy between a SO cellular constrained fluid and a free fluid at the same temperature and volume ( S T ,V ) defines a 'communal entropy'. At low density the communal entropy S T ,V is exactly Nk. Thus, we have a check on the reversibility of the SCSO pathway to RCP. If both the equilibrium SCSO state and metastable hard-sphere glass at RCP are indeed identical, the communal entropy must integrate to zero between y = 0 and and y = y RCP in figure 7. This thermodynamic integration is expressed as whereupon an estimate of the residual entropy of RCP relative to fcc can be obtained from the difference between the SCSO and fccSO p(V ) T equations-of-state. The enthalpy change for any system of hard spheres with no attractions is simply dH = pdV ; because dG = dH −T dS, the entropy difference is obtainable from From figure 8 we see that the two equations-of-state are essentially the same at lower densities. They are also close for all pressures below ∼8.5(σ 3 /kT), i.e. near the cusp of the fccSO plot [63]. Figure 8 shows also that for pressures higher than about 15(σ 3 /kT) the residual volume of the SCSO equationof-state is near constant with any further increase in pressure. At all higher pressures there is effectively a cancellation of the two terms in equation (11). Thus, both the fccSO and SCSO systems have the same compressibilities (or heat capacities) at high pressures. The RCP entropy, relative to fcc at same state point, is obtained by matching numerically the areas as shown in figure 8; when A 3 = A 2 , the residual entropy S 0 = A 1 . The value S 0 = 1.0 ± 0.1 Nk could be sharpened up with more accurate equation-of-state data. Nevertheless, this value agrees with earlier estimates obtained by integrating the heat capacity C p from an equation-of-state for a rapidly quenched glass [40].

Kauzmann's point
It is now more than 30 years since fundamental questions about glass transition phenomenology of hard-sphere fluids were first investigated [40]. At that time, only the virial coefficients up to B 7 were known, and nothing was known about percolation transitions in the hard-sphere fluid. Here, we revisit the hard-sphere metastable fluid region in the light of our new knowledge of higher virial coefficients and the closed-virial equation, the bifurcation at the percolation transition (section 5.2 and figure 5) and a new thermodynamic equation-of-state for the high-density fluid [50]. Using this updated information, we can recalculate the heat capacities of the equilibrium hard-sphere fluid and crystal phases near the freezing transition, and reassess the Kauzmann paradox for hard spheres [69]. This paradox involves the Kauzmann state point, which is defined as the temperature below which a supercooled metastable fluid would have an entropy less than the equilibrium crystal at the same T and p, were it not for the intervention of a glass transition whereupon the heat capacity C p of the supercooled liquid drops suddenly to a value of the crystal.
The HS fcc crystal equation-of-state of Speedy takes free volume form but obeys the self-consistent free volume (SCFV) form [40] where Z = pV /NkT , 'free volume' parameter α o = (1 − V o /V ) and V o is the fcc close-packed volume. Speedy's  (13)), and the crystal equation-of state of Speedy (equation (12)). The high-precision MD pressures of Bannerman [53,54] are shown as open circles. On the top axis for y < y RCP there are 'maximum jammed states' and for y > y RCP a range of hard-sphere 'glassy' states containing frozen crystal-like nucleites.
original values of a, b and c have been sharpened up with the 6-figure accurate MD data of Bannerman [54].
In order to obtain a value for the excess heat capacities of the hard-sphere fluid liquid-like branch for densities above the (available volume) PA percolation threshold we need an equation-of-state. At present, we do not have a highprecision equation comparable to equation (12). For the present purposes, however, an equation-of-state has been found which has a positive pole at the density of Bernal's RCP volume with the form of a remarkably simple free volume equation (no adjustable parameters!): where 'free volume' parameter α B = (1 − V B /V ) and V B is the RCP volume which we obtain from Buffon's constant value (V B = π 2 Nσ 3 /12). The virial expansion equation (4) for the gas-like region, and approximate equation (13) have identical pressures at the percolation transition point y pa = 0.281 [50], and equation (13) reproduces the simulation data in the equilibrium fluid region and metastable fluid region everywhere to within 1% accuracy. Equation (12) gives both the same pressure and heat capacity as the closed virial equation (4) at the accessible volume percolation transition state point The three equations-of-state for hard spheres; gas-like, liquid-like and (fcc) crystal, are shown in the SCFV form in figure 9 for clarity. The equations-of-state (12) and (13) can be used to calculate the heat capacities of the three thermodynamic regions of the hard-sphere model at constant pressure from the exact relationship [40] Figure 10. Excess heat capacities (relative to an ideal gas at the same density and temperature) of the three branches of the hard-sphere fluid and crystal from equations-of-state by equation (14). The solid black line (top left to bottom right) is the virial equation (4); the blue line (near-horizontal) is the crystal phase from the freezing temperature T f to fcc close packing (equation (12)); the red line is the liquid-like region from the percolation temperature T PA to RCP (equation (13)). T K is the Kauzmann point, predicted using the extrapolated virial equation.
where Z is the volume derivative of Z, i.e. (dZ/dV ) p,T . The results are shown in figure 10.
The first observation that we make from figure 10 is that when the Kauzmann point is recalculated, now using 6figure precision crystal data and closed virial equations-ofstate, a result obtained 30 years ago is confirmed unchanged. The resultant, zero entropy difference, Kauzmann point is NkT (K) /pσ 3 = 0.0236 and the corresponding Kauzmann density is y (K) = 0.6362. This reaffirms the coincidence of the RCP density and Kauzmann density noted previously in [40]. Since we do not at present know what if any is the thermodynamic status of the virial equation of state for densities above y PA , this is yet another interesting coincidence to explain.
We can see from the Kauzmann plot in figure 10 that at the freezing transition there is so little difference between the heat capacities of the equilibrium fluid and the crystal that the entropy of fusion, approximately Nk, becomes the residual entropy of the RCP close-packed state relative to the fcc crystal at the same pressure and temperature. This longstanding result is consistent with the residual entropy obtained from figure 8. Thus, we are drawn to a revised conclusion that could not be foreseen 30 years ago [40]: unlike most real liquids, the equilibrium and metastable hard-sphere fluid, as approximately described by equation (13), is continuous in all its derivatives from PA (the available volume percolation transition) to RCP, and shows no Kauzmann catastrophe point.
Nonetheless there remains a question mark about the thermodynamic significance, if any, of the virial equation for densities above the percolation threshold; it is still a very accurate representation for densities up to and, to a small extent, beyond the freezing transition. We know from the form and constants that it contains information about nucleate clusters, and indeed also maximum crystal close packing. To try to resolve this conundrum we look again at ideas that date back 70 years to the Russian scientist Frenkel [70] which suggest that all real systems, especially in the vicinity of phase transitions, will not obey the equations-of-state of the idealized thermodynamic states as set out by J W Gibbs [71]. Rather, the experimental thermodynamic properties can be a hybrid combination of two state functions of competing phases.

Heterophase fluctuations
Frenkel [70] assumed that, in contradiction to the classical thermodynamic treatment of Gibbs, even in the region of thermodynamic stability of a phase A (here the fluid between percolation and freezing), and within the vicinity of phase transitions, the experimental phase is heterogeneous due to the presence of heterophase fluctuations of small nucleates of phase B (here 'solid-like', possibly microcrystalline). These heterophase fluctuations will show up in both simulations and real experiments on finite real systems with surfaces, but not in the general thermodynamic theory of infinite pure, homogeneous single phases in the idealized 'thermodynamic limit' of Gibbs.
These heterophase fluctuations can contribute to anomalous increases in heat capacity in the vicinity of first-order phase transitions: these are well-known experimental and computational phenomena. Both the crystal equation (12) and the virial equation (4) do show abrupt increases in C p at the point of the phase transition for hard spheres. This behaviour may be explained in the case of the virial expansion: we know that it contains information about all clusters in the general configurational integral, but would not show up in the partition function of a thermodynamically ideal pure fluid phase by virtue of the 'maximum term principle' of statistical thermodynamics. The thermodynamic measured average over the probabilities of different phases is equal to the most probable in the limit of an infinite system. Hence we might expect a second bifurcation of the metastable thermodynamic equation-of-state and the virial equation in the vicinity of the freezing transition. Such an anomaly has indeed been suggested by Speedy [72] and Kolafa [73] and observed both in recent simulations [54], and by van Megen and coworkers in real experiments on colloidal hard-sphere suspensions [74][75][76].

Square-well liquids
We have so far focussed entirely on the hard-sphere system, and reviewed compelling evidence that Bernal's RCP is wellcharacterized, well-defined, has a thermodynamic status, and is an appropriate starting point for a theory of the hardsphere 'liquid' branch. In these final three sections, we move on from the hard-sphere system to explore the published evidence relating to its possible significance in understanding the structures of real liquids-i.e. we address how the situation might be changed when we consider the effects of intermolecular attractive forces. Remembering Rowlinson's comment quoted earlier [44] that the attractive forces have 'little direct effect on the structure of a simple liquid', we will not be surprised to conclude that the evidence points to there being indeed a central role for RCP in understanding the structures of real liquids.
When we consider attractive forces, however, any fully successful model of liquids must also be able eventually to explain criticality and critical condensation. In this context it has long since been conjectured [46][47][48][49][50] that the hard-sphere available-volume percolation transition might also be related to the origin of criticality and vapour-liquid coexistence when an attractive perturbation is added. Recent related developments that we discuss below concerning criticality do not conform with the long-held van der Waals hypothesis of continuity of liquid and gas. As these developments offer new insights into longstanding controversies regarding criticality, we proceed to discuss them below within the context of RCP-related concepts originally proposed by Bernal over 50 years ago.

Liquid-vapour coexistence
It follows that if the RCP state of hard spheres has a thermodynamic status that connects it to the hard-sphere fluid at densities above percolation and through the freezing transition in the metastable region, the same must also apply to squarewell liquids. In this section, we see that the evidence that this is indeed the case is already out there in the computer simulation literature of investigations into liquid-gas coexistence.
The square-well potential is a model pair potential between two spheres that can attract each other through a constant potential well depth within a fixed range. For two molecules i and j with hard-core diameter σ , the square-well interaction between them at a separation r is φ ij (r) = 0 for r ij > λσ; −ε for σ > r ij > λσ; ∞ for r ij < σ, where ε is the depth of the attractive well and λσ is its range.
The thermodynamic liquid-vapour coexistence data and critical points of square-well model fluids were first obtained by Vega et al [77] from Gibbs ensemble MC calculations for several values of the square-well width (λ) in the range λ = 1.25σ to 2.0σ . They found that square-well liquids in the coexistence range of temperatures obey the law of rectilinear diameters (LRDs) 5 [78]. They fitted the density curves of the coexisting phases to the following equations which imply obedience to LRD: where ρ c is a 'critical density' in the van der Waals hypothesis and T c the critical temperature. From simulation studies over a range of λ, Vega et al [77] reported results for the critical density for the three constants C 2 , B 0 and β. Setting T = 0 K, the constant mean density in the LRD (ρ m ) is simply where B o now becomes the density of the coexisting 'liquid phase' at absolute zero. It has been observed [63] that the mean value of B o from the literature data averaged over all λ values for which there is data available is 0.60 ± 0.04, and the average value of ρ m = ρ c + C 2 is 0.35 ± 0.1. The within-error equalities of the RCP packing density with the mean value of B 0 , and of half the RCP density (0.318) with the LRD mean density constant ρ m , are thus consistent with the RCP state of the hard-sphere fluid as the metastable amorphous ground state, not only of the hard-sphere fluid, but also of square-well liquids. To echo Rowlinson again [44]: the existence of the attractive forces indeed does have little direct effect on the essential structure of a simple liquid.

The critical region
Though the van der Waals theory of a critical point has been accepted for well over a century as the description of liquidgas criticality, it is worth noting that van der Waals himself recognized in his Nobel lecture [2] that the concept of a critical volume was a weak part in his theory. It should therefore perhaps come as no surprise that an increasing amount of recent published work [63,79,80] presents evidence that the van der Waals hypothetical critical point does not exist on the Gibbs density surface. Building on the extensive calculations of Vega et al [77], a revised phase diagram, consistent with these recent findings, is illustrated in figure 11 for the square-well model fluid (λ = 2). To within the simulation uncertainties, the limiting amorphous ground state of squarewell liquids is the RCP state of hard-spheres as illustrated in figure 5. It can also be observed that, insofar as the simulations permit, the liquid coexistence lines terminate at limiting coexisting maximum liquid density around the region of a percolation transition density. The density of the available volume percolation transition appears to be increased slightly above the unperturbed hard-sphere fluid value (y pa = 0.281) by the attractive-well perturbation.
In addition to one available-volume percolation transition PA, square-well fluids have two excluded-volume percolation transitions. The first is that defined by the hard-sphere excluded volume at the exclusion range 2σ . In addition, the square-well width λ defines a cluster if two or more spheres are linked together within this distance. This leads us to the definition of the 'bonded cluster' percolation transition (PB); it occurs in square-well fluids at a density at which a cluster of bonded atoms first becomes of the same order as system size, i.e. becomes liquid-like. The locus of the bonded cluster percolation transition on the Gibbs density surface has been designated PB. For all square-well fluids with λ < 2σ , which relate to real atomic or molecular fluids, PB occurs at a higher pressure and density than the excluded volume percolation. Figure 11. A phase diagram of a square-well fluid: this example (λ = 2) is plotted from parameterized computer simulation data [77] to within the uncertainties quoted by Vega et al, retaining their critical exponent (β) values but using the present known value of RCP packing fraction i.e. y RCP = 0.6366, and approximate percolation loci [63]. The thin dashed lines are the critical and triple point isobars; solid lines are the gas and liquid coexistence curves. Liquid (blue) bound PA and gas (red) bound PB are the 'available volume' and 'bonded cluster' percolation loci that bound the supercritical liquid and gas phases respectively (thick dashed lines). This transition will first intersect with the PA pressure p pa (T ), and trigger the transition to two coexisting phases of gas and liquid at the critical temperature T c . For the special case of a square-well fluid defined by λ = 2σ , the hard sphere excluded volume and bonded cluster percolation loci become the same. The calculations indicate (see figure 5) that the maximum coexisting gas density coincides with the bonded cluster percolation transition PB at around a packing fraction y pb = 0.05 [49,50]. Figures 5 and 11 indicate that there is no 'continuity of liquid and gas', as hypothesized by van der Waals [1], and parameterized in his legendary equation (1). Rather, a liquid phase spans all temperatures from the metastable amorphous ground state, to supercritical temperatures where it is bounded from a supercritical mesophase by the loci of a percolation transition of the available volume (PA) [78]. Another 'bonded cluster' percolation transition (PB) [63,78] bounds the gas phase. A revision of the properties of liquid argon has revealed a connection between the percolation loci and the metastable spinodals (lines defining a stability limit of a phase). Both PA and PB continue into the two-phase region to define the limits of existence of the metastable liquid and gas phases, respectively. In effect, the spinodals are the extensions of the percolation loci into the liquid and vapour coexistence region [80].
It is evident, therefore, that the RCP state is not just a starting point for the hard-sphere fluid at high density, but also for all model liquids that can be treated theoretically as perturbations of the hard-sphere fluid [4]. The RCP density state bounds the liquid phase diagram, from metastable T * = 0 to equilibrium liquid and beyond T * c all the way to T * = ∞, if crystallization is excluded. Here we have only presented the special case of λ = 2, but the experiment and simulation results [78,79] and theoretical considerations [80] suggests that the Bernal RCP state is the amorphous ground state for all square-well fluids from sticky spheres (λ → σ ) to the limiting mean-field model (λ → ∞).
The existence of a limiting thermodynamic amorphous ground-state cannot be confined to simple model Hamiltonians of hard-sphere and square-well fluids. Next, we consider the relevance to real liquids.

Liquid argon
In section 9 we brought together evidence to show that model square-well liquids have an amorphous ground state at 0 K defined by the extrapolated equation-of-state of the hardsphere fluid. It follows therefore, that if the liquid of a model Hamiltonian has an amorphous ground state, then so also will any real liquid that can be described by a classical Hamiltonian, for example liquid argon. In order to test whether the present experimental data on liquid argon in the equilibrium region is consistent with an extrapolation into the metastable region and a LRD constant consistent with an RCP ground state (see section 9), we can convert it provided we know what the effective hard-sphere diameter of amorphous argon at 0 K is. Clearly, we have no theoretical prescription at present for determining this effective hard-sphere diameter. It seems reasonable to assume, however, that it must lie between the distances of zero attractive potential (σ ) and of the point of zero force (r o ), in accord with the various perturbation theories [45]. Note also that for the Lennard-Jones pair potential r o = 2 (1/6) σ .
The experimental equation-of-state data for liquid argon along the coexistence line is shown in figure 12. Liquid and gas densities, together with coexistence pressures, were reported with 6-figure precision by Gilgen et al [83]. The Gilgen experimental data is here re-plotted in a reduced RCP form. In order to show that the experimental density data of liquid argon is not inconsistent with the idea of an amorphous RCP ground state, a reduced density is calculated as follows. We take the density constant in the LRD obtained by Gilgen et al in SI units of kg m −3 ; defining the mean density at a given T , p as ρ m = (ρ l + ρ v )/2 they obtained for the LRD giving T c = 150.687 K and ρ c = 535.6 kg m −3 . They found that equation (18) corresponds to a straight line through the critical density and the density ρ tr at the triple point temperature (84 K). Setting T = 0 we obtain for the limiting value of ρ m (at T = 0 K) 1.735 × 535.6 = 929.27 kg m −3 . Given the atomic weight and molar volume of atomic argon, we equate this limiting density of Gilgen et al with the RCP packing fraction to see if the resultant effective hard-sphere diameter one obtains is sensibly in the right region of the interatomic pair potential. We obtain a value for this effective argon-argon hard-sphere diameter of σ (eff) = 35.21 nm. This is close to the

The nature of the liquid mesophase
Next, we turn to experimental measurements of the p-V -T properties of liquid argon ( figure 13) not directly in coexistence with the vapour [83]. First we observe that the supercritical region can be divided into three, rather than the classical two, regions as illustrated by the red trapezoid. A gas phase at low density is characterized by a negative rigidity derivative with (d 2 p/dρ 2 ) T < 0; the liquid state by contrast is characterized by a positive rigidity derivative with (d 2 p/dρ 2 ) T > 0. In between the gas and liquid is the proposed mesophase discussed at the end of the previous section 9; we envisage this as a colloidal mixture of gas and liquid with a linear combination rule as a consequence of which (d 2 p/dρ 2 ) T = 0.
A second observation from the experimental isotherms in figure 13 is that it appears to be inconsistent with the coexistence data in figure 12 in the vicinity of the triple point temperature (within the circle at the bottom right of figure 13). This is likely to be a consequence of heterophase fluctuations, in which different experimental conditions can yield different results for the same thermodynamic equilibrium property in the vicinity of phase transitions (see also the discussion on Frenkel's ideas in section 8.2). The isotherms in the vicinity of coexistence of the liquid at low temperature have a slightly higher density than would be expected by the LRDs. Thus we appear to have here an example of the effect of heterophase fluctuations, i.e. 'nucleate clusters' of higher density, which would frustrate any attempt to extrapolate the  [83]; the trapezoid shows the region occupied by the proposed supercritical mesophase. The law of rectilinear diameters (red arrow) is heading for a density of RCP/2 but the coexistence liquid density takes a pre-freezing twist near the triple point (within the red circle). The proposed mesophase is bounded by (on the gas side) the bonded cluster percolation transition PB and (on the liquid side) the available volume percolation transition PA. Reprinted from [83] with permission from Elsevier.
LRD in the region of the freezing temperature to an amorphous ground state. It appears that there is a bifurcation of the experimental liquid y(p, T ) equation-of-state causing a prefreezing compressibility and also heat capacity anomaly that will vitiate the extrapolation. Indeed this non-equilibrium nonanalytic equation of state might be expected to extrapolate not to an amorphous ground state or RCP, but to a crystalline density fcc or hcp.

Bernal's hypercritical line
In his Bakerian lecture to the Royal Society of London in 1962 [7], Bernal made some quite remarkable observations that appear to have anticipated by half a century these developing ideas of criticality discussed above and the proposed supercritical mesophase. Jones and Walker [84] first reported measurements of the specific heat of supercritical fluid argon above the critical temperature. They observed a line of maximum values along an extension of the liquid/vapour line (see figure 14(a) from [84]), and commented on 'how sharply the extrapolated vapourization curve divides the fluid part of the field into liquid-like and vapour-like regions'. In looking for a structural interpretation of this mysterious new line, Bernal commented that it implied there remains a sharp distinction between liquid and gas well above the critical point, and that the critical point itself is only one point on this line which he named the 'hypercritical line'.
We can now re-examine Bernal's hypercritical line and see how it might relate to the concept of the supercritical mesophase bounded by the percolation loci PA and PB, which have been superimposed on the 'parallel' (p, T ) plot of figure 14(b). The critical point occurs where the PA and PB loci coincide while the hypercritical line continues into the mesophase region between PA and PB at higher temperatures and pressures.
With quite remarkable insight, Bernal described graphically ( figure 15) the change in structure of a fluid in going from 'liquid to gas' in the supercritical region  [84]. In his Bakerian lecture [7] Bernal named the line of maximum heat capacity (C p ) that extends beyond the critical point, as the 'hypercritical line'. (b) p-T projection of the liquid-gas coexistence line from measurements of Gilgen et al [81] on which has been superimposed the percolation loci PA and PB [80]. The intersection of PA and PB defines the critical 'divide' at T c , p c whereupon two phases coexist.  [7] showing his supercritical ideas of the transition of molecular arrangements passing from 'expanded' liquid to 'associated gas' through the hypercritical line (dashed). Reproduced from [7] with permission of the Royal Society.
when crossing the hypercritical line at constant temperature. The simple observation he made was that coherence in a liquid was sensitive to volume and not temperature, which is now a well-known phenomenon, whereas a gas has no coherence. His simple explanatory diagram (figure 15) shows not two regions of liquid and gas, but four, with a distinction between coherent 'close packed' and 'expanded' liquid, and the incoherent 'free' and 'associated' gas. The hypercritical line is crossed as the system changes from being 'expanded' liquid to 'associated gas'.
We do not at this stage have detailed structures for the supercritical mesophases of real macroscopic liquids, such as argon, but there can be little doubt that the mesophase is of a colloidal nature [80]. This implies that it must divide into two regions. On the gas side, bounded by the bonded cluster percolation transition loci PB, the colloidal nature will be a dispersed phase of liquid 'microdroplets' and a continuous phase of gas, i.e. a mist. On the liquid side the dispersed phase and continuous phase will be reversed: gas bubbles in continuous liquid, i.e. a foam. At a density (or pressure) intermediate between the percolation loci PA and PB there must therefore be a colloidal inversion from mist to foam. This would involve a greater degree of fluctuations in energy and density, which being exactly related to the heat capacity C p , would result in maximum value loci as seen in figure 14(a). Thus, Bernal's 'hypercritical line' extending the boiling point [7] is also consistent with the new thermodynamic description of criticality that is emerging from recent investigations [79,80,85].

Conclusions
The first part of this topical review explored the light that has been thrown by the recent published literature on the significance for liquid structure of Bernal's inherently simple concept of RCP. Focussing initially on the hardsphere fluid, the evidence points strongly to RCP being not only well-characterized and well-defined, but also having the thermodynamic status of a metastable amorphous ground state.
Whereas the gas phase is described by one equationof-state, i.e. a virial expansion, the liquid state will be described by another that is not in any way connected to the gas equation, and which could be a free volume or cell-cluster expansion. The liquid phase spans the whole pressure and temperature range from an amorphous ground state at zero K, to the supercritical region and beyond.
The evidence we review here suggests that the starting point for the equation-of-state of the hard-sphere fluid in the liquid-like density range is not Mayer's low-density gas virial expansion, nor the fcc crystal, but Bernal's state of RCP.
The relevance of RCP is not, however, limited to the hard sphere or square-well liquids. We can infer that if a model fluid has a reversible path to a metastable amorphous ground state, which is RCP in the case of the square-well model, then real liquids should behave likewise. It seems that there will exist an amorphous ground state approaching absolute zero that terminates a metastable liquid at low temperatures, and that this ground state is defined by extrapolation of the liquid coexistence line to low temperatures, i.e. below the triple point.
It should be possible to obtain an experimental estimate of the residual properties of the amorphous ground state of any liquid from the constants in the LRDs in the equilibrium coexistence region. For simple liquids, whence the linear law is an accurate representation, the limiting density on the metastable coexistence line is 2y m , where y m is the density constant in the LRDs. LRD has a long history [78] in the determination of van der Waals hypothetical 'critical density', which the recent work we have discussed argues does not exist as such [79,80]. The significance of the limiting density constant, however, in its relation to an amorphous ground state with a thermodynamic status, has, until recently, been overlooked [85].
Although the van der Waals picture is widely accepted as an approximate description of the behaviour of real fluids, many have questioned its applicability in the critical region. For example, a review by Levelt-Sengers [86], states 'in a narrow range around the critical point, Van der Waals' equation is fundamentally inadequate'. Moreover it reveals that many early scientists expressed alternative descriptions or reservations of the van der Waals 'continuity of liquid and gas' hypothesis, including such distinguished names as Faraday, Mendeleev, Ramsay, Gouy, Mathias and Cailletet. Moreover, notwithstanding extensive work on the theory of liquid-gas continuity and critical point on the Gibbs density surface [87,88], van der Waals' 'critical point' does not obey the Gibbs phase rule and has no thermodynamic definition. These are serious defects for a concept that is argued to be fundamental to our understanding of liquids, and for 140 years, its existence has remained an unsubstantiated hypothesis.
Looking in detail at the results of extensive calculations on model square-well fluids, radical alternative interpretations have been proposed that suggest a resolution of the conceptual problem of liquid-gas continuity above a 'critical point' of van der Waals. From these ideas emerges a description of the liquid state that is not consistent with the van der Waals idea of continuity of liquid and gas. In contrast to that picture, the reinterpretation discussed above argues rather that the gas and liquid phases in the supercritical region are bounded by percolation transition loci and are well-separated by a colloidal mesophase [63,79,80]. The existence and structure of this mesophase was envisaged by Bernal some 50 years or so ago [7], albeit without knowledge of the thermodynamic implications of percolation phenomena.
The revised phase diagram of simple fluids [79,80,85], based upon percolation transition loci that are determined both from simulations and from experimental density-pressure isotherms of argon in the literature, shows three fluid phases. The supercritical mesophase of a colloidal nature occurs between two percolation transition loci, PA and PB, loci that bound the existences of the liquid and gas phase, respectively. When two percolation transitions have the same pressure, i.e. on intersection of the loci in the p-T plane (e.g. see figure 14(b)), an equilibrium dividing line between the supercritical mesophase and two-phase coexistence is thermodynamically defined. Percolation loci, for T < T c , extend into the gas-liquid coexistence region to become the spinodal lines that limit the existence of metastable gas and liquid phases. These newly proposed liquid-state limits [79,80] now require an equation-of-state for the high-density liquid and supercooled liquid phase that is independent of the gas phase.
Regarding future objectives, we still need to find a theory for the RCP state based upon the statistical geometry that Bernal repeatedly called for over half a century ago; such a theory would likely resolve the puzzling coincidence of the RCP packing fraction with Buffon's constant. There also needs to be a more formal explanation of the RCP residual entropy (relative to fcc) which is close to the ideal gas constant. Further application of high-performance computing is also needed to determine accurately the entire equations-of-state and phase diagrams of all square-well fluids, not least to consolidate these recent findings. Perhaps experiment will be able to verify, and explore the detailed microscopic structures of the 'mist' and 'foam' regions of the mesophase, and the transition between them across the hypercritical line. Moreover, it has not passed our notice that these ideas are also potentially relevant to a huge class of colloidal systems, some of which have important biological relevance: an example might be the coagulation of proteins that are modelled by sticky spheres with a square-well potential [89].
There are various experimental obstacles preventing the determination of the existence and nature of an amorphous ground state. For the likes of liquid argon the problem is spontaneous nucleation and crystallization. In more complex liquids, on supercooling to temperatures well below freezing there is the intervention of a glass transition accompanied by a sudden drop in heat capacity. In these cases, the residual entropies and volumes can be quite small, and it can be difficult to distinguish various degrees of order in predominantly amorphous assemblies. We have also seen that the law of rectilinear diameters as a tool for extrapolating to estimate the density of an amorphous ground state, can be misleading as it fails when there are heterophase fluctuations, as highlighted in the case of argon, for example, in figure 15. There may even be non-analytic behaviour of thermodynamic state functions in the equilibrium liquid close to freezing, as fluctuations that may involve incipient nucleation begin to appear.
This review of the recent science of liquid-state boundaries, whilst resolving queries about the inadequacy of van der Waals hypothesis, also raises more questions and challenges. What now seems clear, however, is that the two original concepts of Bernal, the RCP state of hard spheres, and the 'hypercritical line', will play a significant future role in the theory of liquids.