Modeling the evolution of the melt front under gravity in the presence of a horizontal thermal gradient is a challenging issue, hitherto tackled exclusively with the concepts and tools of computational continuum thermomechanics, too phenomenologically driven to have satisfactory predictive capabilities. Here, we show that this complex phenomenon is amenable to treatment by the methods and tools of Non-Equilibrium Molecular Dynamics (NEMD). To do so, we addressed all the difficulties caused by the necessity of applying suitable boundary conditions and minimizing surface effects so that the bulk behavior of the system in non-equilibrium conditions can be detected. Sufficient adiabatic separation of the time scales permits us to use macroscopically relatively short-but microscopically long enough-time averages to get the macroscopic bulk behavior of the system accurately. To get an adequate signal-to-noise ratio, we had to use an unphysically large value of the gravity. However, we know from NEMD simulations in transport studies that the phenomena produced are stable over many orders of magnitude. In conclusion, our work proves that molecular simulation can be a good tool to study this family of non-equilibrium phenomena, although further work is needed to achieve quantitative predictive capabilities.

Nucleation and growth of methane clathrate hydrates is an exceptional playground to study crystallisation of multi-component, host-guest crystallites when one of the species forming the crystal, the guest, has a higher concentration in the solid than in the liquid phase. This adds problems related to the transport of the low concentration species, here methane. A key aspect in the modelling of clathrates is the water model employed in the simulation. In previous articles, we compared an all-atom force model, TIP4P/Ewald, with a coarse grain one, which is highly appreciated for its computational efficiency. Here, we perform a complementary analysis considering three all-atoms water models: TIP4P/Ewald, TIP4P/ice and TIP5P. A key difference between these models is that the former predicts a much lower freezing temperature. Intuitively, one expects that to lower freezing temperatures of water correspond to lower water/methane-methane gas-clathrate coexistence ones, which determines the degree of supercooling and the degree of supersaturation. Hence, in the simulation conditions, 250 K (500 atm, and fixed methane molar fraction), one expects computational samples made of TIP4P-ice and TIP5P, with a similar freezing temperature (T-f similar to 273 K), to be more supersaturated with respect to the case of TIP4P-Ew (T-f similar to 245 K), and crystallisation to be faster. Surprisingly, we find that while the nucleation rate is consistent with this prediction, growth rate with TIP4P-ice and TIP5P is much slower than with TIP4P-Ew. The latter was attributed to the slower reorientation of water molecules in strong supercooled conditions, resulting in a lower growth rate. This suggests that the freezing temperature is not a suitable parameter to evaluate the adequacy of a water model.[GRAPHICS]

During melting under gravity in the presence of a horizontal thermal gradient, buoyancy-driven convection in the liquid phase affects significantly the evolution of the liquid-solid interface. Due to the obvious engineering interest in understanding and controlling melting processes, fluid dynamicists and applied mathematicians have spent many efforts to model and simulate them numerically. Their endeavors concentrated in the twenty-five years period between the publication of the paper by Brent, Voller & Reid (1988) and that by Mansutti & Bucchignani (2011). The former--and most of the following ones--adopted a phase-field model (where the interface is blurred into a smooth transition zone), while the latter was based on a Stefan-like model with sharp interface. With suitably chosen values of many ad-hoc material and numerical parameters, all of the above simulations were able to attain some agreement with their common benchmark, the melt fronts obtained experimentally by Gau & Viskanta (1986) on a sample of gallium enclosed in a parallelepipedal box with one vertical wall heated. This left unresolved several fine issues, such as whether the elastic response of the solid phase plays a role in determining the shape of the liquid-solid interface.Here, for the first time, we tackle this problem at the atomistic level with a molec- ular dynamics approach. The advantage we gain is that a unique microscopic model describes all of the aggregation states of the molecules, and in particular the solid- liquid interface, without any further assumptions. The price we have to pay is that the hydrodynamical quantities of interest, computed out of the microscopic state using the Irwing & Kirkwood (1950) prescriptions, need to be obtained under gravitational acceleration and thermal gradients much larger than those in real experiments.

Clathrate hydrates are crystalline inclusion compounds wherein a water framework encages small guest atoms/molecules within its cavities. Among the others, methane clathrates are the largest fossil fuel resource still available. They can also be used to safely transport gases and can also form spontaneously under suitable conditions plugging pipelines. Understanding the crystallization mechanism is very important, and given the impossibility of experimentally identifying the atomistic path, simulations played an important role in this field. Given the large computational cost of these simulations, in addition to all-atom force fields, scientists considered coarse-grained water models. Here, we have investigated the effect of coarse-graining, as implemented in the water model mW, on the crystallization characteristics of methane clathrate in comparison with the all-atom TIP4P force field. Our analyses revealed that although the characteristics directly depending on the energetics of the water models are well reproduced, dynamical properties are off by the orders of magnitude. Being crystallization a non-equilibrium process, the altered kinetics of the process results in different characteristics of crystalline nuclei. Both TIP4P and mW water models produce methane clathrate nuclei with some amount of the less stable (in the given thermodynamic conditions) structure II phase and an excess of pentagonal dodecahedral cages over the tetrakaidecahedral ones regarding the ideal ratio in structure I. However, the dependence of this excess on the methane concentration in solution is higher with the former water model, whereas with the latter, the methane concentration in solution dependence is reduced and within the statistical error.

In the framework of the exact factorization of the time-dependent electron-nuclear wave function, we investigate the possibility of solving the nuclear time-dependent Schrödinger equation based on trajectories. The nuclear equation is separated in a Hamilton-Jacobi equation for the phase of the wave function, and a continuity equation for its (squared) modulus. For illustrative adiabatic and nonadiabatic one-dimensional models, we implement a procedure to follow the evolution of the nuclear density along the characteristics of the Hamilton-Jacobi equation. Those characteristics are referred to as quantum trajectories, since they are generated via ordinary differential equations similar to Hamilton's equations, but including the so-called quantum potential, and they can be used to reconstruct exactly the quantum-mechanical nuclear wave function, provided infinite initial conditions are propagated in time.

Derivations of the Jarzynski equality (JE) appear to be quite general, and applicable to any particle system, whether deterministic or stochastic, under equally general perturbations of an initial equilibrium state at given temperatureT. At the same time, the definitions of the quantities appearing in the JE, in particular the work, have been questioned. Answers have been given, but a deeper understanding of the range of phenomena to which the JE applies is necessary, both conceptually and in order to interpret the experiments in which it is used. In fact, domains in which the JE is not applicable have been identified. To clarify the issue, we scrutinize the applicability of the JE to a Hamiltonian particle system in a variable volume. We find that, in this case, the standard interpretation of the terms appearing in the JE is not adequate.

Molecular dynamics (MD) has become one of the most powerful tools of investigation in soft matter. Despite such success, simulations of large molecular environments are mostly run using the approximation of closed systems without the possibility of exchange of matter. Due to the molecular complexity of soft matter systems, an optimal simulation strategy would require the application of concurrent multiscale resolution approaches such that each part of a large system can be considered as an open subsystem at a high resolution embedded in a large coarser reservoir of energy and particles. This paper discusses the current capability and the future perspectives of multiscale adaptive resolution MD methods to satisfy the conceptual principles of open systems and to perform simulations of complex molecular environments in soft matter.

The dihydrogen complex Ru(H2)2H2(P(C5H9)3)2 has been investigated, via ab initio accelerated molecular dynamics, to elucidate the H ligands dynamics and possible reaction paths for H2/H exchange. We have characterized the free energy landscape associated with the H atoms positional exchange around the Ru centre. From the free energy landscape, we have been able to estimate a barrier of 6 kcal mol-1 for the H2/H exchange process. We have also observed a trihydrogen intermediate as a passing state along some of the possible reaction pathways.

Trajectory-based approaches to excited-state, nonadiabatic dynamics are promising simulation techniques to describe the response of complex molecular systems upon photo-excitation. They provide an approximate description of the coupled quantum dynamics of electrons and nuclei trying to access systems of growing complexity. The central question in the design of those approximations is a proper accounting of the coupling electron-nuclei and of the quantum features of the problem. In this paper, we approach the problem in the framework of the exact factorization of the electron-nuclear wavefunction, re-deriving and improving the coupled-trajectory mixed quantum-classical (CT-MQC) algorithm recently developed to solve the exact-factorization equations. In particular, a procedure to include quantum nuclear effects in CT-MQC is derived, and tested on a model system in different regimes.

In various cellular processes, bio-filaments like F-actin and F-tubulin are able to exploit chemical energy associated with polymerization to perform mechanicalwork against an obstacle loaded with an external force. The force-velocity relationship quantitatively summarizes the nature of this process. By a stochastic dynamical model, we give, together with the evolution of a staggered bundle of N-f rigid living filaments facing a loaded wall, the corresponding force-velocity relationship. We compute the evolution of the model in the infinite wall diffusion limit and in supercritical conditions (monomer density reduced by critical density (rho) over cap (1) > 1), and we show that this solution remains valid for moderate non-zero values of the ratio between the wall diffusion and the chemical time scales. We consider two classical protocols: the bundle is opposed either to a constant load or to an optical trap setup, characterized by a harmonic restoring force. The constant load case leads, for each F value, to a stationary velocity V-stat (F; N-f, (rho) over cap (1)) after a relaxation with characteristic time tau(micro)(F). When the bundle (initially taken as an assembly of filament seeds) is subjected to a harmonic restoring force (optical trap load), the bundle elongates and the load increases up to stalling over a characteristic time tau(OT). Extracted from this single experiment, the force-velocity V-OT (F; N-f, (rho) over cap (1)) curve is found to coincide with V-stat (F; N-f, (rho) over cap (1)), except at low loads. We show that this result follows from the adiabatic separation between tau(micro) and tau(OT), i. e., tau(micro) << tau(OT).

A dynamical system submitted to holonomic constraints is Hamiltonian only if considered in the reduced phase space of its generalized coordinates and momenta, which need to be defined ad hoc in each particular case. However, specially in molecular simulations, where the number of degrees of freedom is exceedingly high, the representation in generalized coordinates is completely unsuitable, although conceptually unavoidable, to provide a rigorous description of its evolution and statistical properties. In this paper, we first review the state of the art of the numerical approach that defines the way to conserve exactly the constraint conditions (by an algorithm universally known as SHAKE) and permits integrating the equations of motion directly in the phase space of the natural Cartesian coordinates and momenta of the system. We then discuss in detail SHAKE numerical implementations in the notable cases of Verlet and velocity-Verlet algorithms. After discussing in the same framework how constraints modify the properties of the equilibrium ensemble, we show how, at the price of moving to a dynamical system no more (directly) Hamiltonian, it is possible to provide a direct interpretation of the dynamical system and so derive its Statistical Mechanics both at equilibrium and in non-equilibrium conditions. To achieve that, we generalize the statistical treatment to systems no longer conserving the phase space volume (equivalently, we introduce a non-Euclidean invariant measure in phase space) and derive a generalized Liouville equation describing the ensemble even out of equilibrium. As a result, we can extend the response theory of Kubo (linear and nonlinear) to systems subjected to constraints.

The time reversal invariance of classical dynamics is reconsidered in this paper with specific focus on its consequences for time correlation functions and associated properties such as transport coefficients. We show that, under fairly common assumptions on the interparticle potential, an isolated Hamiltonian system obeys more than one time reversal symmetry and that this entails non trivial consequences. Under an isotropic and homogeneous potential, in particular, eight valid time reversal operations exist. The presence of external fields that reduce the symmetry of space decreases this number, but does not necessarily impair all time reversal symmetries. Thus, analytic predictions of symmetry properties of time correlation functions and, in some cases, even of their null value are still possible. The noteworthy case of a constant external magnetic field, usually assumed to destroy time reversal symmetry, is considered in some detail. We show that, in this case too, some of the new time reversal operations hold, and that this makes it possible to derive relevant properties of correlation functions without the uninteresting inversion of the direction of the magnetic field commonly enforced in the literature.

The paper traces the early stages of Berni Alder's scientific accomplishments, focusing on his contributions to the development of Computational Methods for the study of Statistical Mechanics. Following attempts in the early 50s to implement Monte Carlo methods to study equilibrium properties of many-body systems, Alder developed in collaboration with Tom Wainwright the Molecular Dynamics approach as an alternative tool to Monte Carlo, allowing to extend simulation techniques to non-equilibrium properties. This led to the confirmation of the existence of a phase transition in a system of hard spheres in the late 50s, and was followed by the discovery of the unexpected long-time tail in the correlation function about a decade later. In the late 70s Alder was among the pioneers of the extension of Computer Simulation techniques to Quantum problems. Centered around Alder's own pioneering contributions, the paper covers about thirty years of developments in Molecular Simulation, from the birth of the field to its coming of age as a self-sustained discipline.

We investigate the effects of high solvated-methane concentration on methane-hydrate nucleation at 250 K and 500 atm. We consider solutions at four levels of methane molar fraction in the initial H2O-CH4 solution, ?CH4 = 0.038, 0.044, 0.052, and 0.058, which are higher than (metastable) bulk supersaturation. ?CH4 is controlled independently of the temperature and pressure thanks to the use of special simulation techniques [Phys. Chem. Chem. Phys. 2011, 13, 13177]. These conditions mimic a possible increase of local methane concentration beyond supersaturation induced, for example, by freeze concentration or thermal fluctuations. The nucleation mechanism and kinetics are investigated using the dynamical approach to nonequilibrium molecular dynamics. We demonstrate a hydrate-forming/-ordering process of solvated methane and water molecules in a manner consistent with both the "blob" hypothesis and "cage adsorption hypothesis": the system initially forms an amorphous nucleus at high methane concentration, which then gets ordered, forming the clathrate crystallite. We evaluate nucleation rates using both the methods of the mean first-passage time, i.e., the curve of the average time the system takes to reach a crystalline nucleus of given size, and survival probability, i.e., probability that up to a given time the system has not nucleated yet. We found a dependence of the nucleation rate on initial methane concentration of a form consistent with the dependence of classical nucleation theory rate on supersaturation and determined the relevant parameters of this relation. We found a very rapid increase of nucleation rate with solvated-methane concentration, proving that methane molar fraction, even beyond bulk supersaturation, is key at triggering the homogeneous nucleation of clathrate. We derive a kinetic equation that allows for estimation of the nucleation rate over a wide range of concentration conditions.

We discuss the problem of partitioning a macroscopic system into a collection of independent subsystems. The partitioning of a system into replica-like subsystems is nowadays a subject of major interest in several fields of theoretical and applied physics. The thermodynamic approach currently favoured by practitioners is based on a phenomenological definition of an interface energy associated with the partition, due to a lack of easily computable expressions for a microscopic (i.e. particle-based) interface energy. In this article, we outline a general approach to derive sharp and computable bounds for the interface free energy in terms of microscopic statistical quantities. We discuss potential applications in nanothermodynamics and outline possible future directions.

The time-reversal properties of charged systems in a constant external magnetic field are reconsidered in this paper. We show that the evolution equations of the system are invariant under a new symmetry operation that implies a new signature property for time-correlation functions under time reversal. We then show how these findings can be combined with a previously identified symmetry to determine, for example, null components of the correlation functions of velocities and currents and of the associated transport coefficients. These theoretical predictions are illustrated by molecular dynamics simulations of superionic AgI.

We report a coarse-grained molecular dynamics simulation study of a bundle of parallel actin filaments under supercritical conditions pressing against a loaded mobile wall using a particle-based approach where each particle represents an actin unit. The filaments are grafted to a fixed wall at one end and are reactive at the other end, where they can perform single monomer (de) polymerization steps and push on a mobile obstacle. We simulate a reactive grand canonical ensemble in a box of fixed transverse area A, with a fixed number of grafted filaments Nf, at temperature T and monomer chemical potential m 1. For a single filament case (Nf = 1) and for a bundle of Nf = 8 filaments, we analyze the structural and dynamical properties at equilibrium where the external load compensates the average force exerted by the bundle. The dynamics of the bundle-moving-wall unit are characteristic of an over-damped Brownian oscillator in agreement with recent in vitro experiments by an optical trap setup. We analyze the influence of the pressing wall on the kinetic rates of (de) polymerization events for the filaments. Both static and dynamic results compare reasonably well with recent theoretical treatments of the same system. Thus, we consider the proposed model as a good tool to investigate the properties of a bundle of living filaments.

We establish the statistical mechanics framework for a bundle of N-f living and uncrosslinked actin filaments in a supercritical solution of free monomers pressing against a mobile wall. The filaments are anchored normally to a fixed planar surface at one of their ends and, because of their limited flexibility, they grow almost parallel to each other. Their growing ends hit a moving obstacle, depicted as a second planar wall, parallel to the previous one and subjected to a harmonic compressive force. The force constant is denoted as the trap strength while the distance between the two walls as the trap length to make contact with the experimental optical trap apparatus. For an ideal solution of reactive filaments and free monomers at fixed free monomer chemical potential mu(1), we obtain the general expression for the grand potential from which we derive averages and distributions of relevant physical quantities, namely, the obstacle position, the bundle polymerization force, and the number of filaments in direct contact with the wall. The grafted living filaments are modeled as discrete Wormlike chains, with F-actin persistence length l(p), subject to discrete contour length variations +/- d (the monomer size) to model single monomer (de) polymerization steps. Rigid filaments (l(p) = infinity), either isolated or in bundles, all provide average values of the stalling force in agreement with Hill's predictions F-s(H) = N(f)k(B)T ln(rho(1)/rho(1c))/d, independent of the average trap length. Here rho(1) is the density of free monomers in the solution and rho(1c) its critical value at which the filament does not grow nor shrink in the absence of external forces. Flexible filaments (l(p) < infinity) instead, for values of the trap strength suitable to prevent their lateral escape, provide an average bundle force and an average trap length slightly larger than the corresponding rigid cases (few percents). Still the stalling force remains nearly independent on the average trap length, but results from the product of two strongly L-dependent contributions: the fraction of touching filaments proportional to (< L >(O.T.))(2) and the single filament buckling force proportional to (< L >(O.T.))

Unlike for systems in equilibrium, a straightforward definition of a metastable set in the non-stationary, non-equilibrium case may only be given case-by-case-and therefore it is not directly useful any more, in particular in cases where the slowest relaxation time scales are comparable to the time scales at which the external field driving the system varies. We generalize the concept of metastability by relying on the theory of coherent sets. A pair of sets A and B is called coherent with respect to the time interval [t(1),t(2)] if (a) most of the trajectories starting in A at t(1) end up in B at t(2) and (b) most of the trajectories arriving in B at t(2) actually started from A at t(1). Based on this definition, we can show how to compute coherent sets and then derive finite-time non-stationary Markov state models. We illustrate this concept and its main differences to equilibrium Markov state modeling on simple, one-dimensional examples.

The establishment of thermal diffusion in an Ar-Kr Lennard-Jones mixture is investigated via dynamical non equilibrium molecular dynamics [G. Ciccotti, G. Jacucci, Phys. Rev. Lett. 35, 789 (1975)]. We observe, in particular, the evolution of the density and temperature fields of the system following the onset of the thermal gradient. In stationary conditions, we also compute the Soret coefficient of the mixture. This study confirms that dynamical non equilibrium molecular dynamics is an effective tool to gather information on transient phenomena, even though the full evolution of the mass and energy fluxes associated to the temperature and density fields requires, in this case, a more substantial numerical effort than the one employed here.