Volume 137, Issue 3, 21 July 2012

The systems of openended carbon nanotubes(CNTs) immersed in methanolwater solution are studied by molecular dynamics simulations. For the (6,6) CNT, nearly pure methanol is found to preferentially occupy interior space of the CNT. Even when the mass fraction (MF) of methanol in bulk solution is as low as 1%, the methanol MF within the CNT is still more than 90%. For CNTs with larger diameters, the methanol concentrations within CNTs are also much higher than those outside CNTs. The methanol selectivity decreases with increasing CNT diameter, but not monotonically. From microscopic structural analyses, we find that the primary reason for the high selectivity of methanol by CNTs lies on high preference of methanol in the first solvation shell near the inner wall of CNT, which stems from a synergy effect of the van der Waals interaction between CNT and the methyl groups of methanol, together with the hydrogen bondinginteraction among the liquid molecules. This synergy effect may be of general significance and extended to other systems, such as ethanol aqueous solution and methanol/ethanol mixture. The selective adsorption of methanol over water in CNTs may find applications in separation of water and methanol, detection of methanol, and preservation of methanol purity in fuel cells.
 COMMUNICATIONS


Communication: On the origin of the nonArrhenius behavior in water reorientation dynamics
View Description Hide DescriptionWe combine molecular dynamics simulations and analytic modeling to determine the origin of the nonArrhenius temperature dependence of liquid water's reorientation and hydrogenbond dynamics between 235 K and 350 K. We present a quantitative model connecting hydrogenbond exchange dynamics to local structural fluctuations, measured by the asphericity of Voronoi cells associated with each water molecule. For a fixed local structure the regular Arrhenius behavior is recovered, and the global anomalous temperature dependence is demonstrated to essentially result from a continuous shift in the unimodal structure distribution upon cooling. The nonArrhenius behavior can thus be explained without invoking an equilibrium between distinct structures. In addition, the large width of the homogeneous structural distribution is shown to cause a growing dynamical heterogeneity and a nonexponential relaxation at low temperature.

Communication: Shock adiabat of atomic nitrogen at megabar pressures
View Description Hide DescriptionWe use spherical cellular method combined with a selfconsistent density functional approach (quasizone method) to calculate the band structure and bulk properties of atomic nitrogen at megabar pressures and densities 3.2÷3.6 g/cm^{3}. Thermodynamic functions and shock adiabat calculated by this method correspond to recent measurements, showing a sharp increase in pressure along the shock adiabat in this range of densities.
 Top

 ARTICLES

 Theoretical Methods and Algorithms

Transitiondensityfragment interaction combined with transfer integral approach for excitationenergy transfer via chargetransfer states
View Description Hide DescriptionA transitiondensityfragment interaction (TDFI) combined with a transfer integral (TI) method is proposed. The TDFI method was previously developed for describing electronic Coulomb interaction, which was applied to excitationenergy transfer (EET) [K. J. Fujimoto and S. Hayashi, J. Am. Chem. Soc. 131, 14152 (2009)] and excitoncoupled circular dichroism spectra [K. J. Fujimoto, J. Chem. Phys. 133, 124101 (2010)]. In the present study, the TDFI method is extended to the exchange interaction, and hence it is combined with the TI method for applying to the EET via chargetransfer(CT) states. In this scheme, the overlap correction is also taken into account. To check the TDFITI accuracy, several test calculations are performed to an ethylene dimer. As a result, the TDFITI method gives a much improved description of the electronic coupling, compared with the previous TDFI method. Based on the successful description of the electronic coupling, the decomposition analysis is also performed with the TDFITI method. The present analysis clearly shows a large contribution from the Coulomb interaction in most of the cases, and a significant influence of the CT states at the small separation. In addition, the exchange interaction is found to be small in this system. The present approach is useful for analyzing and understanding the mechanism of EET.

Scalable properties of metal clusters: A comparative study of modern exchangecorrelation functionals
View Description Hide DescriptionThe performance of eight generalized gradient approximation exchangecorrelation (xc) functionals is assessed by a series of scalar relativistic allelectron calculations on octahedral palladium model clusters Pd_{ n } with n = 13, 19, 38, 55, 79, 147 and the analogous clusters Au _{ n } (for n up through 79). For these model systems, we determined the cohesive energies and average bond lengths of the optimized octahedral structures. We extrapolate these values to the bulk limits and compare with the corresponding experimental values. While the wellestablished functionals BP, PBE, and PW91 are the most accurate at predicting energies, the more recent forms PBEsol, VMTsol, and VT{84}sol significantly improve the accuracy of geometries. The observed trends are largely similar for both Pd and Au. In the same spirit, we also studied the scalability of the ionization potentials and electron affinities of the Pd clusters, and extrapolated those quantities to estimates of the work function. Overall, the xc functionals can be classified into four distinct groups according to the accuracy of the computed parameters. These results allow a judicious selection of xc approximations for treating transition metal clusters.

On the accuracy of explicitly correlated coupledcluster interaction energies — have orbital results been beaten yet?
View Description Hide DescriptionThe basis set convergence of weak interaction energies for dimers of noble gases helium through krypton is studied for six variants of the explicitly correlated, frozen geminal coupledcluster singles, doubles, and noniterative triples [CCSD(T)F12] approach: the CCSD(T)F12a, CCSD(T)F12b, and CCSD(T)(F12*) methods with scaled and unscaled triples. These dimers were chosen because CCSD(T) completebasisset (CBS) limit benchmarks are available for them to a particularly high precision. The dependence of interaction energies on the auxiliary basis sets has been investigated and it was found that the default resolutionofidentity sets ccpVXZ/JKFIT are far from adequate in this case. Overall, employing the explicitly correlated approach clearly speeds up the basis set convergence of CCSD(T) interaction energies, however, quite surprisingly, the improvement is not as large as the one achieved by a simple addition of bond functions to the orbital basis set. Bond functions substantially improve the CCSD(T)F12 interaction energies as well. For small and moderate bases with bond functions, the accuracy delivered by the CCSD(T)F12 approach cannot be matched by conventional CCSD(T). However, the latter method in the largest available bases still delivers the CBS limit to a better precision than CCSD(T)F12 in the largest bases available for that approach. Our calculations suggest that the primary reason for the limited accuracy of the largebasis CCSD(T)F12 treatment are the approximations made at the CCSDF12 level and the nonexplicitly correlated treatment of triples. In contrast, the explicitly correlated secondorder MøllerPlesset perturbation theory (MP2F12) approach is able to pinpoint the completebasisset limit MP2 interaction energies of rare gas dimers to a better precision than conventional MP2. Finally, we report and analyze an unexpected failure of the CCSD(T)F12 method to deliver the corecore and corevalence correlation corrections to interaction energies consistently and accurately.

The role of the magnetic orbitals in the calculation of the magnetic coupling constants from multireference perturbation theory methods
View Description Hide DescriptionThe use of multireference perturbation theory (MRPT) for the calculation of the magnetic coupling in binuclear complexes has shown to give poor results if applied on a minimal active space complete active space selfconsistent field (CASSCF) wavefunction. In this work, we identify the origin of this problem in the starting CASSCF orbitals, which are exceedingly localized on the metal atoms. Focusing on the case of antiferromagnetic systems, it is shown that the form of the active orbitals has a dramatic effect on the relative description of the neutral and ionic structures. Finally, a simple and computational inexpensive strategy is proposed for the calculation of a set of magnetic orbitals describing in a more balanced way the neutral and ionic structures. The use of these orbitals, instead the CASSCF ones, in minimal active space MRPT2 calculations leads to a marked improvement of the J values, which become in reasonable agreement with those obtained with the expensive high level difference dedicated configuration interaction approach and with the experimental values.

An algorithm to find minimum freeenergy paths using umbrella integration
View Description Hide DescriptionThe calculation of freeenergy barriers by umbrella sampling and many other methods is hampered by the necessity for an a priori choice of the reaction coordinate along which to sample. We avoid this problem by providing a method to search for saddle points on the freeenergy surface in many coordinates. The necessary gradients and Hessians of the free energy are obtained by multidimensional umbrella integration. We construct the minimum freeenergy path by following the gradient down to minima on the freeenergy surface. The change of free energy along the path is obtained by integrating out all coordinates orthogonal to the path. While we expect the method to be applicable to large systems, we test it on the alanine dipeptide in vacuum. The minima, transition states, and freeenergy barriers agree well with those obtained previously with other methods.

Automatic identification of model reductions for discrete stochastic simulation
View Description Hide DescriptionMultiple time scales in cellular chemical reaction systems present a challenge for the efficiency of stochastic simulation. Numerous model reductions have been proposed to accelerate the simulation of chemically reacting systems by exploiting time scale separation. However, these are often identified and deployed manually, requiring expert knowledge. This is timeconsuming, prone to error, and opportunities for model reduction may be missed, particularly for large models. We propose an automatic model analysis algorithm using an adaptively weighted Petri net to dynamically identify opportunities for model reductions for both the stochastic simulation algorithm and tauleaping simulation, with no requirement of expert knowledge input. Results are presented to demonstrate the utility and effectiveness of this approach.

Treating molecules in arbitrary spin states using the parametric twoelectron reduceddensitymatrix method
View Description Hide DescriptionMinimizing the electronic energy with respect to a parameterized twoelectron reduced density matrix (2RDM) is known as a parametric variational 2RDM method. The parametric 2RDM method with the M 2RDM parametrization [D. A. Mazziotti, Phys. Rev. Lett.101, 253002 (2008)]10.1103/PhysRevLett.101.253002 is extended to treat molecules in arbitrary spin states. Like its singlet counterpart, the M parametric 2RDM method for arbitrary spin states is derived using approximate Nrepresentability conditions, which allow it to capture more correlation energy than coupled cluster with single and double excitations at a lower computational cost. We present energies, optimized bond lengths, potential energy curves, and occupation numbers for a set of molecules in a variety of spin states using the M and K parametric 2RDM methods as well as several wavefunction methods. We show that the M parametric 2RDM method can describe bond breaking of openshell molecules like triplet and singlet and triplet even in the presence of strong correlation. Finally, the computed 2RDMs are shown to be nearly Nrepresentable at both equilibrium and nonequilibrium geometries.

Polaronic discontinuities induced by offdiagonal coupling
View Description Hide DescriptionIn this paper, we study a form of the Holstein molecular crystal model in which the influence of lattice vibrations on the transfers of electronic excitations between neighboring sites (offdiagonal coupling) is taken into account. Using the Toyozawa Ansatz and the Lanczos algorithm, the Holstein Hamiltonian with two types of offdiagonal coupling is studied focusing on a number of analyticity issues in the ground state. For finitesized lattices and antisymmetric coupling, a sequence of discontinuities are found in the polaron energy dispersion, the size of the groundstatephonon cloud, and the linearized von Neumann entropy used to quantify the quantum entanglement between the exciton and the phonons in the ground state. Such behavior is accompanied by a shift of the groundstate crystal momentum from zero to nonzero values as the coupling strength is increased. In the thermodynamic limit, all discontinuities associated with antisymmetric coupling vanish except the one corresponding to the initial departure of the groundstate wavevector from the Brillouin zone center. For the case of symmetric offdiagonal coupling, a smooth crossover is found to exist in all parameters regimes.

Exciton transport in thinfilm cyanine dye Jaggregates
View Description Hide DescriptionWe present a theoretical model for the study of exciton dynamics in Jaggregated monolayers of fluorescent dyes. The excitonic evolution is described by a MonteCarlo wave function approach which allows for a unified description of the quantum (ballistic) and classical (diffusive) propagation of an exciton on a lattice in different parameter regimes. The transition between the ballistic and diffusive regime is controlled by static and dynamic disorder. As an example, the model is applied to three cyanine dye Jaggregates: TC, TDBC, and U3. Each of the moleculespecific structure and excitation parameters are estimated using timedependent density functional theory. The excitondiffusion coefficients are calculated and analyzed for different degrees of film disorder and are correlated to the physical properties and the structural arrangement of molecules in the aggregates. Further, excitontransport is anisotropic and dependent on the initial exciton energy. The upperbound estimation of the excitondiffusion length in the TDBC thinfilm Jaggregate is of the order of hundreds of nanometers, which is in good qualitative agreement with the diffusion length estimated from experiments.

Pair diffusion, hydrodynamic interactions, and available volume in dense fluids
View Description Hide DescriptionWe calculate the pair diffusion coefficient D(r) as a function of the distance r between two hard sphere particles in a dense monodisperse fluid. The distancedependent pair diffusion coefficient describes the hydrodynamic interactions between particles in a fluid that are central to theories of polymer and colloid dynamics. We determine D(r) from the propagators (Green's functions) of particle pairs obtained from molecular dynamics simulations. At distances exceeding ∼3 molecular diameters, the calculated pair diffusion coefficients are in excellent agreement with predictions from exact macroscopic hydrodynamictheory for large Brownian particles suspended in a solvent bath, as well as the Oseen approximation. However, the asymptotic 1/r distance dependence of D(r) associated with hydrodynamic effects emerges only after the pair distance dynamics has been followed for relatively long times, indicating nonnegligible memory effects in the pair diffusion at short times. Deviations of the calculated D(r) from the hydrodynamic models at short distances r reflect the underlying manybody fluid structure, and are found to be correlated to differences in the local available volume. The procedure used here to determine the pair diffusion coefficients can also be used for singleparticle diffusion in confinement with spherical symmetry.

Free energy decomposition analysis of bonding and nonbonding interactions in solution
View Description Hide DescriptionA free energy decomposition analysis algorithm for bonding and nonbonding interactions in various solvated environments, named energy decomposition analysis–polarizable continuum model (EDAPCM), is implemented based on the localized molecular orbitalenergy decomposition analysis (LMOEDA) method, which is recently developed for interaction analysis in gas phase [P. F. Su and H. Li, J. Chem. Phys.130, 074109 (2009)]10.1063/1.3077917. For single determinant wave functions, the EDAPCM method divides the interaction energy into electrostatic, exchange, repulsion, polarization, desolvation, and dispersion terms. In the EDAPCM scheme, the homogeneous solvated environment can be treated by the integral equation formulation of PCM (IEFPCM) or conductorlike polarizable continuum model (CPCM) method, while the heterogeneous solvated environment is handled by the HetCPCM method. The EDAPCM is able to obtain physically meaningful interaction analysis in different dielectric environments along the whole potential energy surfaces. Test calculations by MP2 and DFT functionals with homogeneous and heterogeneous solvation, involving hydrogen bonding, vdW interaction, metalligand binding, cationπ, and ionic interaction, show the robustness and adaptability of the EDAPCM method. The computational results stress the importance of solvation effects to the intermolecular interactions in solvated environments.

On solving the master equation in spatially periodic systems
View Description Hide DescriptionWe present a new method for solving the master equation for a system evolving on a spatially periodic network of states. The network contains 2^{ν} images of a “unit cell” of n states, arranged along one direction with periodic boundary conditions at the ends. We analyze the structure of the symmetrized (2^{ν} n) × (2^{ν} n) rate constant matrix for this system and derive a recursive scheme for determining its eigenvalues and eigenvectors, and therefore analytically expressing the timedependent probabilities of all states in the network, based on diagonalizations of n × n matrices formed by consideration of a single unit cell. We apply our new method to the problem of lowtemperature, lowoccupancy diffusion of xenon in the zeolite silicalite1 using the states, interstate transitions, and transition state theorybased rate constants previously derived by June et al. [J. Phys. Chem. 95 , 8866 (1991)]. The new method yields a diffusion tensor for this system which differs by less than 3% from the values derived previously via kinetic Monte Carlo(KMC) simulations and confirmed by new KMC simulations conducted in the present work. The computational requirements of the new method are compared against those of KMC,numerical solution of the master equation by the Euler method, and direct molecular dynamics. In the problem of diffusion of xenon in silicalite1, the new method is shown to be faster than these alternative methods by factors of about 3.177 × 10^{4}, 4.237 × 10^{3}, and 1.75 × 10^{7}, respectively. The computational savings and ease of setting up calculations afforded by the new method of master equationsolution by recursive reduction of dimensionality in diagonalizing the rate constant matrix make it attractive as a means of predicting longtime dynamical phenomena in spatially periodic systems from atomiclevel information.

Entangled trajectory molecular dynamics in multidimensional systems: Twodimensional quantum tunneling through the Eckart barrier
View Description Hide DescriptionIn this paper, we extend the entangled trajectory molecular dynamics (ETMD) method to multidimensional systems. The integrodifferential form of the evolution equation for the Wigner function is employed, allowing general potentials not represented as a polynomial to be treated. As the example, the method is applied to a twodimensional model of scattering from an Eckart barrier. The results of ETMD are in good agreement with quantum hydrodynamics and exact quantum simulations. By comparing the quantum and classical trajectory in phase space, the quantum tunneling phenomenon is interpreted vividly.

Efficient method to include nuclear quantum effects in the determination of phase boundaries
View Description Hide DescriptionWe developed a methodology to assess nuclear quantum effects in phase boundaries calculations that is based on the dynamical integration of ClausiusClapeyron equation using path integral simulations. The technique employs nonequilibrium simulations that are very efficient. The approach was applied to the calculation of the melting line of Ne in an interval of pressures ranging from 1 to 3366 bar. Our results show a very good agreement with both experimental findings and results from previous calculations. The methodology can be applied to solid and liquid phases, without limitations regarding anharmonicities. The method allows the computation of coexistence lines for wide intervals of pressure and temperature using, in principle, a single simulation.

Scalar fundamental measure theory for hard spheres in three dimensions: Application to hydrophobic solvation
View Description Hide DescriptionHardsphere mixtures provide one a solvable reference system that can be used to improve the density functional theory of realistic molecular fluids. We show how the Kierlik–Rosinberg's scalar version of the fundamental measure density functional theory of hard spheres [E. Kierlik and M. L. Rosinberg, Phys. Rev. A42, 3382 (1990)10.1103/PhysRevA.42.3382], which presents computational advantages with respect to the original Rosenfeld's vectorial formulation or its extensions, can be implemented and minimized in three dimensions to describe fluid mixtures in complex environments. This implementation is used as a basis for defining a molecular density functional theory of water around molecular hydrophobicsolutes of arbitrary shape.

Optimizing conical intersections of solvated molecules: The combined spinflip density functional theory/effective fragment potential method
View Description Hide DescriptionSolvent effects on a potential energy surface crossing are investigated by optimizing a conical intersection (CI) in solution. To this end, the analytic energy gradient has been derived and implemented for the collinear spinflip density functional theory (SFDFT) combined with the effective fragment potential (EFP) solvent model. The new method is applied to the azomethanewater cluster and the chromophore of green fluorescent protein in aqueous solution. These applications illustrate not only dramatic changes in the CI geometries but also strong stabilization of the CI in a polar solvent. Furthermore, the CI geometries obtained by the hybrid SFDFT/EFP scheme reproduce those by the full SFDFT, indicating that the SFDFT/EFP method is an efficient and promising approach for understanding nonadiabatic processes in solution.
 Advanced Experimental Techniques

On the design of experiments for determining ternary mixture free energies from static light scattering data using a nonlinear partial differential equation
View Description Hide DescriptionWe mathematically design sets of static light scattering experiments to provide for modelindependent measurements of ternary liquid mixing free energies to a desired level of accuracy. A parabolic partial differential equation (PDE), linearized from the full nonlinear PDE[D. Ross, G. Thurston, and C. Lutzer, J. Chem. Phys.129, 064106 (Year: 2008)10.1063/1.2937902], describes how data noise affects the free energies to be inferred. The linearized PDE creates a net of spacelike characteristic curves and orthogonal, timelike curves in the composition triangle, and this net governs diffusion of information coming from light scattering measurements to the free energy. Free energy perturbations induced by a light scattering perturbation diffuse along the characteristic curves and towards their concave sides, with a diffusivity that is proportional to the local characteristic curvature radius. Consequently, static light scattering can determine mixing free energies in regions with convex characteristic curve boundaries, given suitable boundary data. The dielectric coefficient is a Lyapunov function for the dynamical system whose trajectories are PDE characteristics. Information diffusion is heterogeneous and systemdependent in the composition triangle, since the characteristics depend on molecular interactions and are tangent to liquidliquid phase separation coexistence loci at critical points. We find scaling relations that link free energy accuracy, total measurement time, the number of samples, and the interpolation method, and identify the key quantitative tradeoffs between devoting time to measuring more samples, or fewer samples more accurately. For each total measurement time there are optimal sample numbers beyond which more will not improve free energy accuracy. We estimate the degree to which manypoint interpolation and optimized measurement concentrations can improve accuracy and save time. For a modest light scattering setup, a sample calculation shows that less than two minutes of measurement time is, in principle, sufficient to determine the dimensionless mixing free energy of a nonassociating ternary mixture to within an integrated error norm of 0.003. These findings establish a quantitative framework for designing light scattering experiments to determine the Gibbs free energy of ternary liquid mixtures.

Mathematical and computational aspects of quaternary liquid mixing free energy measurement using light scattering
View Description Hide DescriptionWe provide a mathematical and computational analysis of light scatteringmeasurement of mixing free energies of quaternary isotropic liquids. In previous work, we analyzed mathematical and experimental design considerations for the ternary mixture case [D. Ross, G. Thurston, and C. Lutzer, J. Chem. Phys.129, 064106 (2008)10.1063/1.2937902;C. Wahle, D. Ross, and G. Thurston, J. Chem. Phys.137, 034201 (2012)10.1063/1.4731694]. Here, we review and introduce dimensionfree general formulations of the fully nonlinear partial differential equation(PDE) and its linearization, a basis for applying the method to composition spaces of any dimension, in principle. With numerical analysis of the PDE as applied to the light scattering implied by a test free energy and dielectric gradient combination, we show that values of the Rayleigh ratio within the quaternary composition tetrahedron can be used to correctly reconstruct the composition dependence of the free energy. We then extend the analysis to the case of a finite number of data points, measured with noise. In this context the linearized PDE describes the relevant diffusion of information from light scattering noise to the free energy. The fully nonlinear PDE creates a special set of curves in the composition tetrahedron, collections of which form characteristics of the nonlinear and linear PDEs, and we show that the information diffusion has a timelike direction along the positive normals to these curves. With use of Monte Carlo simulations of light scattering experiments, we find that for a modest laboratory light scattering setup, about 100–200 samples and 100 s of measurement time are enough to be able to measure the mixing free energy over the entire quaternary composition tetrahedron, to within an error norm of 10^{−3}. The present method can help quantify thermodynamics of quaternary isotropic liquid mixtures.