Percolation of polyatomic species with the presence of impurities

In this paper, the percolation of (a) linear segments of size k and (b) k-mers of different structures and forms deposited on a square lattice contaminated with previously adsorbed impurities have been studied. The contaminated or diluted lattice is built by randomly selecting a fraction of the elements of the lattice (either bonds or sites) which are considered forbidden for deposition. Results are obtained by extensive use of finite size scaling theory. Thus, in order to test the universality of the phase transition occurring in the system, the numerical values of the critical exponents were determined. The characteristic parameters of the percolation problem are dependent not only on the form and structure of the k-mers but also on the properties of the lattice where they are deposited. A phase diagram separating a percolating from a nonpercolating region is determined as a function of the parameters of the problem. A comparison between random site and random bond percolation in the presence of impurities on the lattice is presented.     

Critical behavior of long straight rigid rods on two-dimensional lattices: theory and Monte Carlo simulations

The critical behavior of long straight rigid rods of length k (k-mers) on square and triangular lattices at intermediate density has been studied. A nematic phase, characterized by a big domain of parallel k-mers, was found. This ordered phase is separated from the isotropic state by a continuous transition occurring at an intermediate density theta(c). Two analytical techniques were combined with Monte Carlo simulations to predict the dependence of theta(c) on k, being theta(c)(k) proportional to k(-1). The first involves simple geometrical arguments, while the second is based on entropy considerations. Our analysis allowed us also to determine the minimum value of k (k(min) = 7), which allows the formation of a nematic phase on a triangular lattice.     

Monte Carlo simulations of two-dimensional hard core lattice gases

Monte Carlo simulations are used to study lattice gases of particles with extended hard cores on a two-dimensional square lattice. Exclusions of one and up to five nearest neighbors (NN) are considered. These can be mapped onto hard squares of varying side length, lambda (in lattice units), tilted by some angle with respect to the original lattice. In agreement with earlier studies, the 1NN exclusion undergoes a continuous order-disorder transition in the Ising universality class. Surprisingly, we find that the lattice gas with exclusions of up to second nearest neighbors (2NN) also undergoes a continuous phase transition in the Ising universality class, while the Landau-Lifshitz theory predicts that this transition should be in the universality class of the XY model with cubic anisotropy. The lattice gas of 3NN exclusions is found to undergo a discontinuous order-disorder transition, in agreement with the earlier transfer matrix calculations and the Landau-Lifshitz theory. On the other hand, the gas of 4NN exclusions once again exhibits a continuous phase transition in the Ising universality class-contradicting the predictions of the Landau-Lifshitz theory. Finally, the lattice gas of 5NN exclusions is found to undergo a discontinuous phase transition.     

Cluster algorithm for nonadditive hard-core mixtures

In this paper, we present a cluster algorithm for the numerical simulations of nonadditive hard-core mixtures. This algorithm allows one to simulate and equilibrate systems with a number of particles two orders of magnitude larger than previous simulations. The phase separation for symmetric binary mixtures is studied for different nonadditivities as well as for the Widom-Rowlinson model [B. Widom and J. S. Rowlinson, J. Chem. Phys. 52, 1670 (1970)] in two and three dimensions. The critical densities are determined from finite size scaling. The critical exponents for all the nonadditivities are consistent with the Ising universality class.     

Liquid-liquid phase transitions in supercooled water studied by computer simulations of various water models

Liquid-liquid and liquid-vapor coexistence regions of various water models were determined by Monte Carlo (MC) simulations of isotherms of density fluctuation-restricted systems and by Gibbs ensemble MC simulations. All studied water models show multiple liquid-liquid phase transitions in the supercooled region: we observe two transitions of the TIP4P, TIP5P, and SPCE models and three transitions of the ST2 model. The location of these phase transitions with respect to the liquid-vapor coexistence curve and the glass temperature is highly sensitive to the water model and its implementation. We suggest that the apparent thermodynamic singularity of real liquid water in the supercooled region at about 228 K is caused by an approach to the spinodal of the first (lowest density) liquid-liquid phase transition. The well-known density maximum of liquid water at 277 K is related to the second liquid-liquid phase transition, which is located at positive pressures with a critical point close to the maximum. A possible order parameter and the universality class of liquid-liquid phase transitions in one-component fluids are discussed.     

On the definition of a Monte Carlo model for binary crystal growth

We show that consistency of the transition probabilities in a lattice Monte Carlo (MC) model for binary crystal growth with the thermodynamic properties of a system does not guarantee the MC simulations near equilibrium to be in agreement with the thermodynamic equilibrium phase diagram for that system. The deviations remain small for systems with small bond energies, but they can increase significantly for systems with large melting entropy, typical for molecular systems. These deviations are attributed to the surface kinetics, which is responsible for a metastable zone below the liquidus line where no growth occurs, even in the absence of a 2D nucleation barrier. Here we propose an extension of the MC model that introduces a freedom of choice in the transition probabilities while staying within the thermodynamic constraints. This freedom can be used to eliminate the discrepancy between the MC simulations and the thermodynamic equilibrium phase diagram. Agreement is achieved for that choice of the transition probabilities yielding the fastest decrease of the free energy (i.e., largest growth rate) of the system at a temperature slightly below the equilibrium temperature. An analytical model is developed, which reproduces quite well the MC results, enabling a straightforward determination of the optimal set of transition probabilities. Application of both the MC and analytical model to conditions well away from equilibrium, giving rise to kinetic phase diagrams, shows that the effect of kinetics on segregation is even stronger than that predicted by previous models.     

The condensation and ordering of models of empty liquids

We consider a simple model consisting of particles with four bonding sites ("patches"), two of type A and two of type B, on the square lattice, and investigate its global phase behavior by simulations and theory. We set the interaction between B patches to zero and calculate the phase diagram as the ratio between the AB and the AA interactions, ε(AB)*, varies. In line with previous work, on three-dimensional off-lattice models, we show that the liquid-vapor phase diagram exhibits a re-entrant or "pinched" shape for the same range of ε(AB)*, suggesting that the ratio of the energy scales--and the corresponding empty fluid regime--is independent of the dimensionality of the system and of the lattice structure. In addition, the model exhibits an order-disorder transition that is ferromagnetic in the re-entrant regime. The use of low-dimensional lattice models allows the simulation of sufficiently large systems to establish the nature of the liquid-vapor critical points and to describe the structure of the liquid phase in the empty fluid regime, where the size of the "voids" increases as the temperature decreases. We have found that the liquid-vapor critical point is in the 2D Ising universality class, with a scaling region that decreases rapidly as the temperature decreases. The results of simulations and theoretical analysis suggest that the line of order-disorder transitions intersects the condensation line at a multi-critical point at zero temperature and density, for patchy particle models with a re-entrant, empty fluid, regime.     

Energetics driving the short-range order in CuxPd1-x/Ru(0001) monolayer surface alloys

The energetics determining the distinct short-range order in two-dimensional (2D) monolayer Cu(x)Pd(1-x) surface alloys on a Ru(0001) substrate were investigated by Monte Carlo simulations and density functional theory calculations. Using a 2D lattice gas Hamiltonian based on effective pair interaction (EPI) parameters, the EPIs were derived for different Cu concentrations with Monte Carlo (MC) simulations by comparing with the atomic distributions obtained from atomic resolution STM images and the related Warren-Cowley short-range order parameters (Hoster et al., Phys. Rev. B, 2006, 73 165413). The ground state structures and mixing energies at 0 K derived from these EPIs agree well with mixing energies determined from DFT calculations of different ordered surface alloys. Additional MC simulations yield rather low transition temperatures which explain the absence of ordered 2D phases in the experiments. The consequences of our findings for the use of alloy surfaces and surface alloys as model systems for adsorption and catalytic reaction studies are discussed.     

Prediction of a structural transition in the hard disk fluid

Starting from the second equilibrium equation in the BBGKY hierarchy under the Kirkwood superposition closure, we implement a new method for studying the asymptotic decay of correlations in the hard disk fluid in the high density regime. From our analysis and complementary numerical studies, we find that exponentially damped oscillations can occur only up to a packing fraction η(?)～0.718, a value that is in substantial agreement with the packing fraction, η～0.723, believed to characterize the transition from the ordered solid phase to a dense fluid phase, as inferred from Mak's Monte Carlo simulations [Phys. Rev. E 73, 065104 (2006)]. Next, we show that the same method of analysis predicts that the exponential damping of oscillations in the hard sphere fluid becomes impossible when λ=4nπσ(3)[1+H(1)]≥34.81, where H(1) is the contact value of the correlation function, n is the number density, and σ is the sphere diameter in exact agreement with the condition, λ≥34.8, which is first reported in a numerical study of the Kirkwood equation by Kirkwood et al. [J. Chem. Phys. 18, 1040 (1950)]. Finally, we show that our method confirms the absence of any structural transition in hard rods for the entire range of densities below close packing.     

Unimolecular decomposition of formic acid in the gas phase--on the ratio of the competing reaction channels

The thermal decomposition of formic acid was reinvestigated in the gas phase using two types of shock tubes. It was confirmed that the unimolecular decomposition proceeds through a main channel of dehydration (k1) and a minor decarboxylation channel (k2). This result is in good agreement with our previous study (J. Chem. Phys. 1984, 80, 4989). Furthermore, it was confirmed that the dehydration process is in the second-order region and that the decarboxylation is in the falloff region, in the temperature range of 1300-2000 K and over the total density of (0.5-2.5) x 10(-5) mol cm(-3). The experimental ratios between the two channels, k2/k1, are compared with those of theoretical calculations by conventional transition state theory and the Rice-Ramsperger-Kassel-Marcus theory.     

Simulation study of two-dimensional phase transitions of argon on graphite surface and in slit micropores

Molecular simulation has been increasingly used in the analysis and modeling of gas adsorption on open surfaces and in porous materials because greater insight could be gained from such a study. In case of homogeneous surfaces or pore walls the adsorption behavior is often complicated by the order–disorder transition. It is shown in our previous publications (Ustinov and Do, Langmuir 28:9543–9553, 2012a; Ustinov and Do, Adsorption 19:291–304, 2013) that once an ordered molecular layer has been formed on the surface, the lattice constant depends on the simulation box size, which requires adjusting the box dimensions parallel to the surface for each value of loading. It was shown that this can be accomplished with the Gibbs–Duhem equation, which results in decreasing lattice constant with an increase of the amount adsorbed. The same feature is expected to be valid for gas adsorption in narrow pores, but this has not been analyzed in the literature. This study aims at an extension of our approach to adsorption in slit graphitic pores using kinetic Monte Carlo method (Ustinov and Do, J Colloid Interface Sci 366:216–223, 2012b). The emphasis rests on the thermodynamic analysis of the two-dimensional (2D) ordering transition and state of the ordered phase; if the ordered phase exists in narrow slit pores, simulation with constant volume box always leads to erroneous results, for example, seemingly incompressible adsorbed phase. We proposed a new approach that allows for modeling thermodynamically consistent adsorption isotherms, which can be used as a basis for further refinement of the pore size distribution analysis of nanoporous materials.     

Computational predictions of stable 2D arrays of bidisperse particles

We study computationally the stability of various 2D arrays of bidisperse mixtures of stabilized nanoparticles through a melting simulation employing the Metropolis algorithm for determining surface diffusion. In our previous work [Langmuir 2004, 20, 9408], we studied computationally the stability of bidispersed monolayers of thiol-stabilized gold nanoparticles with a size ratio (sigma) of 0.375. We found that interparticle forces were essential to stabilize the LS (the two-dimensional NaCl analogue) lattice at the experimentally determined surface coverage. In this paper, we extend our study to determine the conditions necessary to form stable LS(2), LS(4), and LS(6) lattices, which have yet to be observed. Using a simple design rule that involves matching the distances between either large-large particles and large-small particles or large-small particles and small-small particles to correspond to the respective potential minima leads to predictions for size ratios that will form each desired lattice, given other parameters characterizing the systems' physical properties. We predict and verify computationally LS(2), LS(4), and LS(6) lattices at relatively low surface coverages. Additional simulations show that the LS, LS(2), and LS(6) lattices are indeed stable structures at their predicted surface coverage, whereas the LS(4) lattice is a metastable structure; however, a modest increase in the surface coverage of the LS(4) lattice converts it to a stable rather than long-lived metastable structure. This study may be used as a guide for experimentalists in their search for these novel structures.     

On the existence of a third-order phase transition beyond the Andrews critical point: a molecular dynamics study

The possibility of the existence of a gas-liquid third order phase transition for fluids is becoming a subject of growing interest. Experimental work suggests its existence for specific systems while recent theoretical models claim its universality. In this work, we employ Molecular Dynamics and investigate the third-order phase transition beyond the Andrews critical point by treating a system of Lennard-Jones particles along three isotherms. Two partial derivatives of the Gibbs free energy are measured, namely the molar constant pressure heat capacity and isothermal compressibility. The convergence of these simulations with respect to the system size as well as the cut-off radius is carefully checked. The obtained results show that partial derivatives certainly do not present sharp cusp singularities at the maxima, and actually suggest that there are no singularities at all. On these basis we then conclude that a third-order phase transition in the considered temperature region: T? ≥ 1.36 may indeed not exist.     

Finite-size scaling analysis on the phase transition of a ferromagnetic polymer chain model

The finite-size scaling analysis method is applied to study the phase transition of a self-avoiding walking polymer chain with spatial nearest-neighbor ferromagnetic Ising interaction on the simple cubic lattice. Assuming the scaling M2(T,n) = n(-2beta/nu)[phi0 + phi1n(1/nu)(T-T(c)) + O(n(2/nu)(T-T(c))2)] with the square magnetization M2 as the order parameter and the chain length n as the size, we estimate the second-order phase-transition temperature T(c) = 1.784 J/k(B) and critical exponents 2beta/nu approximately 0.668 and nu approximately 1.0. The self-diffusion constant and the chain dimensions (R2) and (S2) do not obey such a scaling law.     

Clustering in the absence of attractions: density functional theory and computer simulations

We numerically investigate the formation of stable clusters of overlapping particles in certain systems interacting via purely repulsive, bounded pair potentials. In close vicinity of a first-order phase transition between a disordered and an ordered structure, clusters are encountered already in the fluid phase which then freeze into crystals with multiply occupied lattice sites. These hyper-crystals are characterized by a number of remarkable features that are in clear contradiction to our experience with harshly repulsive systems: upon compression, the lattice constant remains invariant, leading to a concomitant linear growth in the cluster population with density; further, the freezing and melting lines are to high accuracy linear in the density-temperature plane, and the conventional indicator that announces freezing, that is, the Hansen-Verlet value of the first peak of the structure factor, attains for these soft systems much higher values than for their hard-matter counterparts. Our investigations are based on the generalized exponential model of index 4 (i.e., Phi(r) approximately exp[-(r/sigma)4]). The properties of the phases involved are calculated via liquid state theory and classical density functional theory. Monte Carlo simulations for selected states confirm the theoretical results for the structural and thermodynamic properties of the system. These numerical data, in turn, fully corroborate an approximate theoretical framework that was recently put forward to explain the clustering phenomenon for systems of this kind (Likos, C. N.; Mladek, B. M.; Gottwald, D.; Kahl, G. J. Chem. Phys. 2007, 126, 224502).     

A molecular dynamics study on universal properties of polymer chains in different solvent qualities. Part I. A review of linear chain properties

This paper investigates the conformational and scaling properties of long linear polymer chains. These investigations are done with the aid of Monte Carlo (MC) and molecular dynamics (MD) simulations. Chain lengths that comprise several orders of magnitude to reduce errors of finite size scaling, including the effect of solvent quality, ranging from the athermal limit over the theta-transition to the collapsed state of chains are investigated. Also the effect of polydispersity on linear chains is included which is an important issue in the real fabrication of polymers. A detailed account of the hybrid MD and MC simulation model and the exploited numerical methods is given. Many results of chain properties in the extrapolated limit of infinite chain lengths are documented and universal properties of the chains within their universality class are given. An example of the difference between scaling exponents observed in actual solvents and those observed in the extremes of "good solvents" and "theta-solvents" in simulations is provided by comparing simulation results with experimental data on low density polyethylene. This paper is concluded with an outlook on the extension of this study to branched chain systems of many different branching types.     

Magnetic tuning of all-organic binary alloys between two stable radicals

Mixtures of 2-(4,5,6,7-tetrafluorobenzimidazol-2-yl)-4,4,5,5-tetramethyl-4,5-dihydro-1H-imidazole-3-oxide-1-oxyl (F4BImNN) and 2-(benzimidazol-2-yl)-4,4,5,5-tetramethyl-4,5-dihydro-1H-imidazole-3-oxide-1-oxyl (BImNN) crystallize as solid solutions (alloys) across a wide range of binary compositions. (F4BImNN)(x)(BImNN)((1-x)) with x < 0.8 gives orthorhombic unit cells, while x ≥ 0.9 gives monoclinic unit cells. In all crystalline samples, the dominant intermolecular packing is controlled by one-dimensional (1D) hydrogen-bonded chains that lead to quasi-1D ferromagnetic behavior. Magnetic analysis over 0.4-300 K indicates ordering with strong 1D ferromagnetic exchange along the chains (J/k = 12-22 K). Interchain exchange is estimated to be 33- to 150-fold weaker, based on antiferromagnetic ordered phase formation below Néel temperatures in the 0.4-1.2 K range for the various compositions. The ordering temperatures of the orthorhombic samples increase linearly as (1 - x) increases from 0.25 to 1.00. The variation is attributed to increased interchain distance corresponding to decreased interchain exchange, when more F4BImNN is added into the orthorhombic lattice. The monoclinic samples are not part of the same trend, due to the different interchain arrangement associated with the phase change.     

Accelerating atomic-level protein simulations by flat-histogram techniques

Flat-histogram techniques provide a powerful approach to the simulation of first-order-like phase transitions and are potentially very useful for protein studies. Here, we test this approach by implicit solvent all-atom Monte Carlo (MC) simulations of peptide aggregation, for a 7-residue fragment (GIIFNEQ) of the Cu/Zn superoxide dismutase 1 protein (SOD1). In simulations with 8 chains, we observe two distinct aggregated/non-aggregated phases. At the midpoint temperature, these phases coexist, separated by a free-energy barrier of height 2.7 k(B)T. We show that this system can be successfully studied by carefully implemented flat-histogram techniques. The frequency of barrier crossing, which is low in conventional canonical simulations, can be increased by turning to a two-step procedure based on the Wang-Landau and multicanonical algorithms.     

Free volume hypothetical scanning molecular dynamics method for the absolute free energy of liquids

The hypothetical scanning (HS) method is a general approach for calculating the absolute entropy, S, and free energy, F, by analyzing Boltzmann samples obtained by Monte Carlo (MC) or molecular dynamics (MD) techniques. With HS applied to a fluid, each configuration i of the sample is reconstructed by gradually placing the molecules in their positions at i using transition probabilities (TPs). With our recent version of HS, called HSMC-EV, each TP is calculated from MC simulations, where the simulated particles are excluded from the volume reconstructed in previous steps. In this paper we remove the excluded volume (EV) restriction, replacing it by a "free volume" (FV) approach. For liquid argon, HSMC-FV leads to an improvement in efficiency over HSMC-EV by a factor of 2-3. Importantly, the FV treatment greatly simplifies the HS implementation for liquids, allowing a much more natural application of the method for MD simulations. Given the success and popularity of MD, the present development of the HSMD method for liquids is an important advancement for HS methodology. Results for the HSMD-FV approach presented here agree well with our HSMC and thermodynamic integration results. The efficiency of HSMD-FV is equivalent to HSMC-EV. The potential use of HSMC(MD)-FV in protein systems with explicit water is discussed.     

Molecular simulation study of cavity-generated instabilities in the superheated Lennard-Jones liquid

Previous equilibrium-based density-functional theory (DFT) analyses of cavity formation in the pure component superheated Lennard-Jones (LJ) liquid [S. Punnathanam and D. S. Corti, J. Chem. Phys. 119, 10224 (2003); M. J. Uline and D. S. Corti, Phys. Rev. Lett. 99, 076102 (2007)] revealed that a thermodynamic limit of stability appears in which no liquidlike density profile can develop for cavity radii greater than some critical size (being a function of temperature and bulk density). The existence of these stability limits was also verified using isothermal-isobaric Monte Carlo (MC) simulations. To test the possible relevance of these limits of stability to a dynamically evolving system, one that may be important for homogeneous bubble nucleation, we perform isothermal-isobaric molecular dynamics (MD) simulations in which cavities of different sizes are placed within the superheated LJ liquid. When the impermeable boundary utilized to generate a cavity is removed, the MD simulations show that the cavity collapses and the overall density of the system remains liquidlike, i.e., the system is stable, when the initial cavity radius is below some certain value. On the other hand, when the initial radius is large enough, the cavity expands and the overall density of the system rapidly decreases toward vaporlike densities, i.e., the system is unstable. Unlike the DFT predictions, however, the transition between stability and instability is not infinitely sharp. The fraction of initial configurations that generate an instability (or a phase separation) increases from zero to unity as the initial cavity radius increases over a relatively narrow range of values, which spans the predicted stability limit obtained from equilibrium MC simulations. The simulation results presented here provide initial evidence that the equilibrium-based stability limits predicted in the previous DFT and MC simulation studies may play some role, yet to be fully determined, in the homogeneous nucleation and growth of embryos within metastable fluids.